SOREL: A Stochastic Algorithm for Spectral Risks Minimization

## Abstract

Abstract The spectral risk has wide applications in machine learning, especially in real-world decision-making, where people are not only concerned with models’ average performance. By assigning different weights to the losses of different sample points, rather than the same weights as in the empirical risk, it allows the model’s performance to lie between the average performance and the worst-case performance. In this paper, we propose SOREL, the first stochastic gradient-based algorithm with convergence guarantees for the spectral risk minimization. Previous algorithms often consider adding a strongly concave function to smooth the spectral risk, thus lacking convergence guarantees for the original spectral risk. We theoretically prove that our algorithm achieves a near-optimal rate of O ~ ⁢ ( 1 / ϵ ) ~ 𝑂 1 italic-ϵ \widetilde{O}(1/\sqrt{\epsilon}) over~ start_ARG italic_O end_ARG ( 1 / square-root start_ARG italic_ϵ end_ARG ) in terms of ϵ italic-ϵ \epsilon italic_ϵ . Experiments on real datasets show that our algorithm outperforms existing algorithms in most cases, both in terms of runtime and sample complexity.

## 1 Introduction

In modern machine learning, model training heavily relies on minimizing the empirical risk. This ensures that the model have high average performance.
Given a training set of $n$ sample points $\mathcal{D}=\{{\bm{x}_{i}}\}_{i=1}^{n}\subset\mathcal{X}$, the goal of the empirical risk minimization is to solve

$$ $R(\bm{w})=(1/n)\sum_{i=1}^{n}\ell_{i}(\bm{w}).$ $$

Here, $\bm{w}\in\mathbb{R}^{d}$ is the parameter vector of the model, $\ell_{i}(\bm{w})=\ell(\bm{w},\bm{x}_{i})$ is the loss of the $i$-th sample, and $\ell:\mathbb{R}^{d}\times\mathcal{X}\to\mathbb{R}$ is the loss function. However, as machine learning models are deployed in real-world scenarios, the evaluation metrics for these models become more diverse, including factors such as fairness or risk aversion.

In this paper, we consider a generalized aggregation loss function: the spectral risk, which is of the form

$$ $R_{\bm{\sigma}}(\bm{w})=\sum_{i=1}^{n}\sigma_{i}\ell_{[i]}(\bm{w}).$ $$

Here $\ell_{[1]}(\cdot)\leq\dots\leq\ell_{[n]}(\cdot)$ denotes the order statistics of the empirical loss distribution, and $0\leq\sigma_{1}\leq\dots\leq\sigma_{n},\sum_{i=1}^{n}\sigma_{i}=1.$

In form, the spectral risk penalizes the occurrence of extreme losses by assigning greater weights to extreme losses. When $\sigma_{i}=1/n$, the spectral risk reduces to the empirical risk. When $\sigma_{n}=1$ and $\sigma_{i}=0$ for $i=1,\dots,n-1$, the spectral risk becomes the maximum loss. Therefore, the spectral risk measures the model’s performance between the average case and the worst case. By assigning different values to $\sigma_{i}$, the spectral risk encompasses a wide range of aggregated loss functions that have broad applications in fields such as machine learning and finance.
Common spectral risks include Conditional Value at Risk (CVaR) or the average of top-k loss [, ], Exponential Spectral Risk Measure (ESRM) [], and Extremal Spectral Risk (Extremile) []. Their specific forms are shown in Table 1 [].

**Table 1: Different spectral risks and the corresponding weights.**
| Spectral Risks | Parameter | $\sigma_{i}$ |
| --- | --- | --- |
| $\alpha$-CVaR | $0<\alpha<1$ | $\begin{cases}\frac{1}{n\alpha},&i>\lceil n(1-\alpha)\rceil\\<br>1-\frac{\lfloor n\alpha\rfloor}{n\alpha},&\lfloor n(1-\alpha)\rfloor<i<\lceil n<br>(1-\alpha)\rceil\\<br>0,&\text{otherwise}\end{cases}$ |
| $\rho$-ESRM | $\rho>0$ | $e^{-\rho}\left(e^{\rho\frac{i}{n}}-e^{\rho\frac{i-1}{n}}\right)/\left(1-e^{-<br>\rho}\right)$ |
| $r$-Extremile | $r\geq 1$ | $\left(\frac{i}{n}\right)^{r}-\left(\frac{i-1}{n}\right)^{r}$ |

Despite the broad applications of spectral risks, optimization methods for spectral risks are still limited. In particular, for large-scale optimization problems, there is currently a lack of stochastic algorithms with convergence guarantees for the spectral risk minimization.
Indeed, the weight of each sample point depend on the entire training set, introducing complex dependencies and thus making the optimization process challenging.
Existing algorithms either abandon the convergence guarantee to the minimum of the spectral risk problem due to the difficulty of obtaining unbiased subgradient estimates [, ], or turn to optimize the smooth regularized spectral risk [, ], which lacks convergence guarantees for the original spectral risk. Given the widespread application of the spectral risk in machine learning and the lack of stochastic algorithms for the spectral risk minimization, we are committed to developing stochastic algorithms with convergence guarantees for the spectral risk minimization.

#### Our Contributions.

In this paper, we study Stochastic Optimization for Spectral Risks with trajectory Stabilization (SOREL). i) We propose SOREL, the first stochastic algorithm with convergence guarantees for the spectral risk minimization. In particular, SOREL stabilizes the trajectory of the primal variable to the optimal point. ii) Theoretically, we prove that SOREL achieves a near-optimal rate of $\widetilde{O}(1/\sqrt{\epsilon})$ in terms of $\epsilon$ for spectral risks with a strongly convex regularizer. This matches the known lower bound of $\Omega(1/\sqrt{\epsilon})$ in the deterministic setting []. iii) Experimentally, SOREL outperforms existing baselines in most cases, both in terms of runtime and sample complexity.

## 2 Related work

#### Statistical Properties of the Spectral Risk.

As a type of risk measures, the spectral risk assigns higher weights to the tail distribution and has been profoundly studied in the financial field [, , ]. Recently, statistical properties of the spectral risk have been investigated by many works in the field of learning theory. Specifically, [, ] demonstrated that the discrete form of spectral risks converges to the spectral risk of the overall distribution, controlled by the Wasserstein distance. [, , ] also considered the learning bounds of spectral risks.

#### Applications.

The spectral risk is widely applied in various fields of finance and machine learning. In some real-world tasks, the worst-case loss is as important as the average-case loss, such as medical imaging [] or robotics []. The spectral risk minimization can be viewed as a risk-averse learning strategy. In the domain of fair machine learning, different subgroups are classified by sensitive features (e.g., gender and race). Subgroups with higher losses may be treated unfairly in decision-making. By penalizing samples with higher losses, the model’s performance across different subgroups can meet certain fairness criteria [], such as demographic parity [] and equalized odds [].
In the field of distributionally robust optimization, the sample distribution may deviate from a uniform distribution, which can be modeled by reweighting the samples []. [] adopts the polytope as the uncertainty set of the shifted distribution, which is similar to the form of the spectral risk minimization.

In practical applications, we can choose different types of spectral risks based on actual needs. For example, CVaR is popular in the context of portfolio optimization [], as well as reinforcement learning [, ]. Other applications of spectral risks include reducing test errors and mitigating the impact of outliers [, , ], to name a few.

#### Existing Optimization Methods.

There have been many algorithms to optimize CVaR, a special case of spectral risks, including both deterministic [] and stochastic algorithms [, ]. For the spectral risk, deterministic methods such as subgradient methods have convergence guarantees, although they are considered algorithms with slow convergence rate. Xiao et al. [] proposed an Alternating Direction Method of Multipliers (ADMM) type method for the minimization of the rank-based loss.
Other deterministic methods reformulate this problem into a minimax problem Thekumparampil et al. [], Hamedani and Aybat [], Khalafi and Boob []. However, these methods require calculating
$O(n)$ function values and gradients at each iteration, posing significant limitations when solving large-scale problems.

Stochastic algorithms for solving the spectral risk minimization problems are still limited. Some algorithms forgo convergence to the true minimum of the spectral risk [, ]. Other methods modify the objective to minimize a smooth approximation of the spectral risk by adding a strongly concave term with a coefficient $\nu$ [, ]. The smaller $\nu$ is, the closer the approximation is to the original spectral risk. [] proposed the Prospect algorithm and proved that it achieves linear convergence for any
$\nu>0$. Furthermore, if the loss of each sample is different at the optimal point, then the optimal value of the smooth approximation of the spectral risk is the same as the optimal value of the original spectral risk
as long as
$\nu$ is below a certain positive threshold. However, in practice, these conditions are difficult to verify. The convergence of these algorithms still lacks guarantees for the original spectral risk minimization.
Other stochastic algorithms include Hamedani and Jalilzadeh [], Yan et al. [], which have a slower convergence rate of $O(1/\epsilon)$ in terms of $\epsilon$. In this paper, we propose SOREL for the original spectral risk minimization problems and achieve a near-optimal convergence rate in terms of
$\epsilon$.

## 3 Algorithm

We consider stochastic optimization of the spectral risk combined with a strongly convex regularizer:

$$ $\displaystyle\min_{\bm{w}}\underbrace{\sum_{i=1}^{n}\sigma_{i}\ell_{[i]}(\bm{w })}_{R_{\bm{\sigma}}(\bm{w})}+g(\bm{w}).$ (1) $$

Firstly, we make the basic assumption about the individual loss function $\ell_{i}$ and the regularizer $g$.

###### Assumption 1 .

The individual loss function $\ell_{i}:\mathbb{R}^{d}\to\mathbb{R}$ is convex, $G$-Lipschitz continuous and $L$-smooth for all $i\in\{1,\dots,n\}$. The regularizer $g:\mathbb{R}^{d}\to\mathbb{R}\cup\{\infty\}$ is proper, lower semicontinuous and $\mu$-strongly convex.

Assumption 1 is common in the literature on stochastic optimization [, ], especially in the field of the spectral risk minimization [, , , , ]. Both the logistic loss and the least-square loss within bounded sublevel sets satisfy this assumption. The $L_{2}$ regularization, which is one of the most common regularizers in machine learning, satisfies the assumption on $g$.

###### Assumption 1 .

### 3.1 Challenges of Stochastic Optimization for Spectral Risks

In this section, we describe the challenges of the spectral risk minimization problem and the techniques to solve them.

#### Biases of Stochastic Subgradient Estimators.

From convex analysis [, Lemma 10], we know that the subgradient of $R_{\bm{\sigma}}$ is

$$ $\partial R_{\bm{\sigma}}(\bm{w})=\operatorname{Conv}\left\{\bigcup_{\pi}\left \{\sum_{i=1}^{n}\sigma_{i}\nabla\ell_{\pi(i)}(\bm{w}):\ell_{\pi(1)}(\bm{w}) \leq\dots\leq\ell_{\pi(n)}(\bm{w})\right\}\right\},$ $$

where $\operatorname{Conv}$ denotes the convex hull of a set, and $\pi$ is a permutation that arranges $\ell_{1},\dots,\ell_{n}$ in ascending order. Note that $R_{\bm{\sigma}}(\bm{w})$ is non-smooth. Indeed, when there exist $i\neq j$ such that $\ell_{i}(\bm{w})=\ell_{j}(\bm{w})$, $\partial R_{\bm{\sigma}}(\bm{w})$ contains multiple elements.

The subgradient of $R_{\bm{\sigma}}$ is related to the ordering of $\ell_{1},...,\ell_{n}$. We cannot obtain an unbiased subgradient estimator of $\partial R_{\bm{\sigma}}$ if we use only a mini-batch with $m$ $(m<n)$ sample points. For example, when $m=1$, we randomly sample $i$ uniformly from $\{1,\dots,n\}$. The subgradient estimator $\nabla\ell_{i}(\bm{w})$ is unbiased only if $\sigma_{i}=1/n$.
For general $\sigma$, unfortunately, to obtain an unbiased subgradient estimator of $\partial R_{\bm{\sigma}}$,
we have to compute $n$ loss function values and then determine the ranking of $\ell_{i}$ among the $n$ losses (or the weight corresponding to the $i$-th sample point). However, computing $O(n)$ losses at each step is computationally heavy. To remedy this, we next design an algorithm that first use a minimax reformulation of problem (1) and then alternates updates $\bm{\lambda}$ and $\bm{w}$ using a primal-dual method.

Equivalently, we can rewrite $R_{\bm{\sigma}}(\bm{w})$ in the following form

$$ $\displaystyle R_{\bm{\sigma}}(\bm{w})=\max_{\bm{\lambda}\in\Pi_{\bm{\sigma}}} \sum_{i=1}^{n}\lambda_{i}\ell_{i}(\bm{w}),$ (2) $$

where $\Pi_{\bm{\sigma}}=\{\Pi\bm{\sigma}:\Pi\mathbf{1}=\mathbf{1},\Pi^{\top}\mathbf{
1}=\mathbf{1},\Pi\in[0,1]^{n\times n}\}$ is the permutahedron associated with $\bm{\sigma}$, i.e., the convex hull of all permutations of $\bm{\sigma}$, and $\mathbf{1}$ is the all-one vector []. Then Problem (1) can be rewritten as

$$ $\displaystyle\min_{\bm{w}}\max_{\bm{\lambda}\in\Pi_{\bm{\sigma}}}L(\bm{w},\bm{ \lambda})=\sum_{i=1}^{n}\lambda_{i}\ell_{i}(\bm{w})+g(\bm{w}).$ (3) $$

Next, we use a primal-dual method to solve Problem (3). Specifically, we iteratively update $\bm{w}$ and $\bm{\lambda}$:

$$ $\bm{\lambda}_{k+1}=\underset{\bm{\lambda}\in\Pi_{\bm{\sigma}}}{\arg\max}\sum_{ i=1}^{n}\lambda_{i}\ell_{i}(\bm{w}_{k})-\frac{1}{2\eta_{k}}\|\bm{\lambda}-\bm{ \lambda}_{k}\|^{2},$ (4) $$

$$ $\bm{w}_{k+1}=\underset{\bm{w}}{\arg\min}P_{k}(\bm{w}):=\sum_{i=1}^{n}\lambda_{ i,k+1}\ell_{i}(\bm{w})+g(\bm{w})+\frac{1}{2\tau_{k}}\|\bm{w}-\bm{w}_{k}\|^{2}.$ (5) $$

Steps (4) and (5) can be seen as alternatingly solving the min problem and the max problem in (3) with proximal terms.

#### Stabilizing the Optimization Trajectory.

To update $\bm{\lambda}_{k+1}$, one may naturally think of solving Problem (2): $\bm{\lambda}_{k+1}=\arg\max_{\bm{\lambda}\in\Pi_{\bm{\sigma}}}\sum_{i=1}^{n}
\lambda_{i}\ell_{i}(\bm{w}_{k})$, similar to methods in Mehta et al. [, ] with smoothing coefficient $\nu=0$. However, since Problem (2) is merely convex, the solution
$\bm{\lambda}$ lacks continuity with respect to $\bm{w}$, that is, a small change in $\bm{w}$ could lead to a large change in $\bm{\lambda}$. Indeed, it is often the case that there are multiple optimal solutions for (2) when there exists $i\neq j$ such that $\ell_{i}(\bm{w})=\ell_{j}(\bm{w})$, and in this case, an arbitrary small perturbation of $\bm{w}$ will lead to a different value of $\lambda_{i}$.
As shown in Figure 1, this can cause
$\bm{w}$ to oscillate near points where some losses are the same and prevents the convergence of the algorithm. We also provide a toy example in Appendix C to further illustrate this difficulty. Therefore, the proximal term $\frac{1}{2\eta_{k}}\|\bm{\lambda}-\bm{\lambda}_{k}\|^{2}$ is added in (4) to prevent excessive changes in $\bm{\lambda}$ and stabilize the trajectory of the primal variable, where $\eta_{k}>0$ controls the extent of its variation.

#### Stochastic Optimization for the Primal Variable.

We use a stochastic algorithm to approximately solve (5). Through the minimax reformulation in (5), we avoid directly calculating the stochastic subgradient of $R_{\bm{\sigma}}(\bm{w})$, which requires computing all loss function values to obtain the corresponding sample weight $\lambda_{i}$. Additionally, since $\bm{\lambda}_{k+1}$ is fixed, the finite sum part of the objective function in (5) is smooth, allowing us to use variance reduction (VR), a commonly used technique in stochastic optimization [, , , ], to accelerate our stochastic algorithm. In contrast, since $R_{\bm{\sigma}}(\bm{w})$ is non-smooth, as previously mentioned, VR cannot be used to directly solve Problem (1). For smooth convex functions in the form of the finite sum, many methods such as SVRG [], SAGA [], and SARAH [] can enable stochastic methods to achieve the convergence rate of deterministic methods. We apply the proximal
stochastic gradient descent with a generalized VR method inspired by SVRG to approximately solve (5), which will be presented in Section 3.2 in detail. Thanks to its strong convexity, Problem (5) can be solved efficiently.

Similar to (4), we add a proximal term $\frac{1}{2\tau_{k}}\|\bm{w}-\bm{w}^{k}\|^{2}$ in (5) where $\tau_{k}>0$ is the proximal parameter. The proximal parameter $\tau_{k}$ is crucial for the convergence proof of our algorithm. By carefully choosing $\tau_{k}=O(1/k)$, the updates of $\bm{w}$ become more stringent as the algorithm progresses, and SOREL can achieve a near optimal rate of $\widetilde{O}(1/\sqrt{\epsilon})$ in terms of $\epsilon$.

Figure: Figure 1: The level set plot of 2D least-square regression with primal-dual optimization trajectories described in Section 3.1. The max subproblem does not have a proximal term (left) or has a proximal term (right). The min subproblem does not have a proximal term. The black star represents the optimal point. The sample points are obtained by projecting the yacht dataset onto $\mathbb{R}^{2}$ using PCA.
Refer to caption: x1.png

### 3.2 The SOREL Algorithm

Our proposed algorithm SOREL is summarized in Algorithm 1. The specific values for the parameters $\theta_{k},\eta_{k},\tau_{k}$ and $m_{k}$ in Algorithm 1 will be given in Section 4. In Line 1 the algorithm initializes $\bm{\lambda}_{0}$ by solving Problem (2). In Lines 7-14, the algorithm computes the stochastic gradient and update $\bm{w}$ for a fixed $\bm{\lambda}$, as described in Section 3.1. Additionally, we compute the full gradient of $\bm{w}$ every $m_{k}$ updates to reduce the variance. In Lines 3-4, we update $\bm{\lambda}$. Note that we replace $\ell_{i}(\bm{w}_{k})$ with $\ell_{i}(\bm{w}_{k})+\theta_{k}\left(\ell_{i}(\bm{w}_{k})-\ell_{i}(\bm{w}_{k-1
})\right)$ to accelerate the algorithm. This can be seen as a momentum term, a widely used technique in smooth optimization [, , , ], where $\theta_{k}>0$ is the momentum parameter.

Figure: Algorithm 1 SOREL

Define the proximal operator $\operatorname{prox}_{h}(\bar{\bm{x}}):=\arg\min_{\bm{x}}h(\bm{x})+\frac{1}{2}
\|\bm{x}-\bar{\bm{x}}\|^{2}$ for a function $h$. In Line 14, we apply the proximal
stochastic gradient descent step. We assume that $\operatorname{prox}_{g+\frac{1}{2}\|\cdot\|^{2}}(\cdot)$ is easy to compute, which is the case for many commonly used regularizers $g$, such as the
$l_{1}$ norm. If $g$ is differentiable, we can replace the proximal stochastic gradient with stochastic gradient: $\bm{w}_{k,t+1}=\bm{w}_{k,t}-\alpha\left(\bm{d}_{k,t}+\frac{1}{\tau_{k}}\left(
\bm{w}_{k,t}-\bm{w}_{k}\right)+\nabla g(\bm{w}_{k,t})\right).$
This will not affect the convergence or convergence rate of the algorithm as long as $\nabla g$ is Lipschitz continuous and the step size $\alpha$ is small enough.
In Line 4, we need to compute the projection onto $\Pi_{\bm{\sigma}}$. For an ordered vector, projecting onto the permutahedron takes $O(n)$ operations using the Pool Adjacent Violators Algorithm (PAVA)[]. In SOREL, we need to first sort $n$ elements of the projected vector and then compute the projection onto $\Pi_{\bm{\sigma}}$, which takes a total of $O(n\log n)$ operations.

In practice, we set $T_{k}$ and $m_{k}$ to $n$ in Lines 7 and 8, meaning the algorithm updates $\bm{\lambda}$ once it traverses the training set.
We also set the reference point $\bar{\bm{w}}$ and the output $\bm{w}_{k+1}$ in Lines 9 and 16 to be the last vector of the previous epoch rather than the average vector, as with most practical algorithms [, , , , ].
In this way, SOREL only requires computing the full batch gradient once for each update of $\bm{\lambda}$, and becomes single-loop in Lines 7- 15. This makes the algorithm more concise and parameters easier to tune.

## 4 Theoretical Analysis

For convenience, we consider that $T_{k}$ (will be determined in Theorem 1) is large enough so that $\bm{w}_{k}$ is a $\delta_{k}$-optimal solution of $P_{k}(\bm{w})$, that is, $\mathbb{E}_{k}P_{k}(\bm{w}_{k+1})-\min_{\bm{w}}P_{k}(\bm{w})\leq\delta_{k}$. Here, $\mathbb{E}_{k}$ represents the conditional expectation with respect to the random sample points used to compute $\bm{w}_{k+1}$ given $\bm{w}_{k},\dots,\bm{w}_{0}$. Then, we can provide a one-step analysis of the outer loop of SOREL. We use $L(\bm{w},\bm{\lambda})=\bm{\lambda}^{\top}\bm{\ell}(\bm{w})+g(\bm{w})$ in the analysis for simplicity.

###### Lemma 1 .

Suppose Assumption 1 holds. Let $\{\bm{w}_{k}\}$ and $\{\bm{\lambda}_{k}\}$ be the sequences generated by Algorithm 1. Then for any $\bm{w}\in\mathbb{R}^{d}$ and $\bm{\lambda}\in\Pi_{\bm{\sigma}}$, the following inequality holds,

$$ $\displaystyle\mathbb{E}_{k}\,\left\{L(\bm{w}_{k+1},\bm{\lambda})-L(\bm{w},\bm{ \lambda}_{k+1})\right\}$ (6) $\displaystyle\leq$ $\displaystyle\mathbb{E}_{k}\left\{\langle\bm{\lambda}-\bm{\lambda}_{k+1},\bm{ \ell}(\bm{w}_{k+1})\rangle+\frac{1}{2\eta_{k}}\left[\|\bm{\lambda}-\bm{\lambda }_{k}\|^{2}-\|\bm{\lambda}-\bm{\lambda}_{k+1}\|^{2}-\|\bm{\lambda}_{k+1}-\bm{ \lambda}_{k}\|^{2}\right]\right.$ $\displaystyle+\frac{1}{2\tau_{k}}\left[\|\bm{w}-\bm{w}_{k}\|^{2}-\|\bm{w}-\bm{ w}_{k+1}\|^{2}-\|\bm{w}_{k+1}-\bm{w}_{k}\|^{2}\right]-\frac{\mu}{2}\|\bm{w}- \bm{w}_{k+1}\|^{2}$ $\displaystyle+\left.\langle\bm{v}_{k},\bm{\lambda}_{k+1}-\bm{\lambda}\rangle+ \delta_{k}+\sqrt{\frac{(\tau_{k}^{-1}+\mu)\delta_{k}}{2}}(\|\bm{w}-\bm{w}_{k+1 }\|^{2}+1)\right\}.$ $$

Next, we try to telescope the terms on the right hand side of (6) by multiplying each term by $\gamma_{k}$.
By choosing appropriate parameters in Algorithm 1 to satisfy some conditions (will be discussed in Appendix A), we can ensure that the adjacent terms indexed by $k=0,\dots,K-1$ can be canceled out during summation. Then we can achieve the desired convergence result.

###### Theorem 1 .

Suppose Assumption 1 holds.
Set $\gamma_{k}=k+1,~{}\eta_{k}=\frac{\mu(k+1)}{8G^{2}},~{}\theta_{k}=\frac{k}{k+1}
,~{}\tau_{k}=\frac{4}{\mu(k+1)},~{}\delta_{k}=\min\left(\frac{\mu}{8(k+5)},\mu
(k+1)^{-6}\right)$,
the step-size $\alpha=\frac{1}{12L}$, $m_{k}=\frac{384L}{(k+5)\mu}+2$ and $T_{k}=O(m_{k}\log\frac{1}{\delta_{k}})$ in Algorithm 1. Then we have

$$ $\mathbb{E}\|\bm{w}_{K}-\bm{w}^{\star}\|^{2}=O\left(\frac{G^{2}}{\mu^{2}K^{2}} \right).$ $$

###### Corollary 1 .

Under the same conditions in Theorem 1, we obtain an output $\bm{w}_{K}$ of Algorithm 1 such that $\mathbb{E}\|\bm{w}_{K}-\bm{w}^{\star}\|^{2}\leq\epsilon$ in a total sample complexity of $O\left(n\frac{G}{\mu\sqrt{\epsilon}}\log\frac{G}{\mu^{2}\sqrt{\epsilon}}+\frac
{L}{\mu}\log\frac{G}{\mu\sqrt{\epsilon}}\log\frac{G}{\mu^{2}\sqrt{\epsilon}}\right)$.

Our algorithm achieves a near-optimal convergence rate of $\widetilde{O}(1/\sqrt{\epsilon})$ in terms of $\epsilon$, which matches the lower bound of $\Omega(1/\sqrt{\epsilon})$ in the deterministic setting up to a logarithmic term []. This is the first near-optimal method for solving the spectral risk minimization. Previously, Mehta et al. [, ] add a strongly concave term with respect to $\bm{\lambda}$ in $L(\bm{w},\bm{\lambda})$ and achieve a linear convergence rate for the perturbed problem.
One may set the coefficient of the strongly concave term $\nu$ to $O(\epsilon)$, obtaining an $\epsilon$-optimal solution for the original spectral risk minimization problem. However, this approach has drawbacks: it leads to a worse sample complexity of $\widetilde{O}(1/\epsilon)$ [] or even $\widetilde{O}(1/\epsilon^{3})$ []; additionally, to achieve an $\epsilon$-optimal solution, the step size would need to be set to $O(\epsilon)$, resulting in very small steps that perform poorly in practice. In contrast, SOREL’s step size is independent of $\epsilon$.

###### Remark 1 .

In Lines 9 and 16 of Algorithm 1, we set the reference point $\bar{\bm{w}}$ and the output $\bm{w}_{k+1}$ to the average of the previous epoch. Instead, we can also set them to be the last vector of the previous epoch, which aligns with practical implementation. For theoretical completeness, we may compute the full gradient $\bar{\bm{g}}$ in Line 10 at each step $t$ with probability $p$ instead of once per epoch (every $m_{k}$ steps), as done in [, , ]. However, these methods are beyond the scope of this paper.

###### Lemma 1 .

###### Theorem 1 .

###### Corollary 1 .

###### Remark 1 .

## 5 Experiments

In this section, we compare our proposed algorithm SOREL with existing baselines for solving the spectral risk minimization problem.

We consider the least squares regression problem using a linear model with $L_{2}$ regularization. We adopt a wide variety of spectral risks, including ESRM, Extremile, and CVaR. Five tabular regression benchmarks are used for the least squares loss: yacht [], energy [], concrete [], kin8nm [], power [].
We compare the suboptimality versus passes (the number of samples divided by $n$) and runtime.
The suboptimality is defined as

$$ $\text{Suboptimality}(\bm{w}_{k})=\frac{R_{\bm{\sigma}}(\bm{w}_{k})+g(\bm{w}_{k })-R_{\bm{\sigma}}(\bm{w}^{\star})-g(\bm{w}^{\star})}{R_{\bm{\sigma}}(\bm{w}_{ 0})+g(\bm{w}_{0})-R_{\bm{\sigma}}(\bm{w}^{\star})-g(\bm{w}^{\star})},$ $$

where $\bm{w}^{\star}$ is calculated by L-BFGS [].

Baseline methods include SGD [, ] with a minibatch size of 64, LSVRG [], and Prospect []. Note that although both LSVRG and Prospect add a strongly concave term with coefficient $\nu$ to smooth the original spectral risk, they have been observed to exhibit linear convergence for the original spectral risk minimization problem in practice without the strongly concave term [, ]. Consequently, we set $\nu=0$ in our experiments. By [], if the losses at the optimal point are different from each other, then as long as $\nu$ is below a certain positive threshold, the optimal solutions of the smoothed and original spectral risks minimization problems are the same. For LSVRG, we set the length of an epoch to $n$. For SOREL, we set $T_{k}=m_{k}=n$, and in all experiments, the values of $\tau_{k}$ are fixed and $\theta_{k}=k/(k+1)$. Thus, our algorithm has only two hyperparameters $\eta_{k}$ and the learning rate $\alpha$ to tune. We apply stochastic gradient descent to solve (5) instead of proximal stochastic gradient descent.
We also empirically compare different algorithms with minibatching, where the batch size for all algorithms are set to 64.
Detailed experimental settings are provided in Appendix B.

Figure: Figure 2: Suboptimality of spectral risks for different algorithms without mini-batching. The $x$-axis represents the effective number of samples used by the algorithm divided by $n$ (odd columns) or CPU time (even columns). Each row corresponds to the same dataset, and each column corresponds to the same type of the spectral risk.
Refer to caption: x2.png

#### Results.

Figure: Figure 3: Suboptimality of spectral risks for different algorithms with mini-batching. The $x$-axis represents the effective number of samples used by the algorithm divided by $n$ (odd columns) or CPU time (even columns).
Refer to caption: x3.png

Figure 2 compares the training curves of our method with other baselines across various datasets and the spectral risk settings.
In terms of sample complexity and runtime, SOREL outperforms other baselines in most cases; SOREL also achieve comparable results in the kin8nm dataset.
In the power dataset, the sample complexity of Prospect is better than that of SOREL. However, the runtime of SOREL is significantly shorter than that of Prospect due to the fact that Prospect needs the calculation of projections onto the permutahedron with $O(n)$ operations each step.
As expected, SGD fails to converge due to its inherent bias [, ]. Although Mehta et al. [] discusses the equivalence of minimizing the smoothed spectral risk and the original spectral risk when losses at the optimal point are different from each other, we find that LSVRG and Prospect often fail to reach the true optimal point, indicating limitations of these methods. In contrast, SOREL converges to the true optimal point in all settings (suboptimality less than 0 means the solution’s accuracy is higher than L-BFGS).

Additionally, in Figure 3, we present the empirical results of the algorithms with minibatching. Minibatching has a significant improvement on the convergence rate of all the algorithms, even though we have not proven the convergence of them.
Similar to what is shown in Figure 2,
SGD, LSVEG and Prospect fail to converge to the true optimal points, especially in the first two datasets.
SOREL converges to the optimal solutions in all settings, and achieve the best or competitive results, in terms of sample complexity, runtime, or both, except in the setting of CVaR and power dataset. Still SOREL performs competitively for an accuracy of $10^{-7}$ in this setting.

## 6 Conclusion

We have proposed SOREL, the first stochastic algorithm with convergence guarantees for the spectral risk minimization problems. We have proved that SOREL achieves a near-optimal rate of $\widetilde{O}(1/\sqrt{\epsilon})$ in terms of $\epsilon$. In experiments, SOREL outperforms existing baselines in terms of sample complexity and runtime in most cases. The SOREL algorithm with minibatching also empirically achieves significant rate improvements.

Future work includes exploring convergence of SOREL for nonconvex problems, and investigating broader applications of the spectral risk in areas such as fairness and distributionally robust optimization.

## Appendix A Proofs

First, we provide an auxiliary lemma. This is an extension of [, Lemma 8].

###### Lemma 2 .

Let $\bar{\bm{x}}$ be an $\epsilon$-approximate solution of $\min_{x}\{g(\bm{x})+\frac{\lambda}{2}\|\bm{x}-\hat{\bm{x}}\|^{2}\}$ in expectation, where $g:\mathbb{R}^{d}\to\mathbb{R}$ is $\mu$-strongly convex, $\mu\geq 0$. Then

$$ $\displaystyle\mathop{\mathbb{E}}\left\{g(\bar{\bm{x}})-g(\bm{x})\right\}\leq$ $\displaystyle\mathbb{E}\left\{\frac{\lambda}{2}\left[\|\bm{x}-\hat{\bm{x}}\|^{ 2}-\|\bm{x}-\bar{\bm{x}}\|^{2}-\|\hat{\bm{x}}-\bar{\bm{x}}\|^{2}\right]-\frac{ \mu}{2}\|\bm{x}-\bar{\bm{x}}\|^{2}\right\}$ $\displaystyle+\sqrt{\frac{(\lambda+\mu)\epsilon}{2}}\mathbb{E}\|\bar{\bm{x}}- \bm{x}\|^{2}+\sqrt{\frac{(\lambda+\mu)\epsilon}{2}}+\epsilon.$ $$

###### Proof.

Let $\bm{x}^{\star}=\arg\min\limits_{x}\{g(\bm{x})+\frac{\lambda}{2}\|\bm{x}-\hat{
\bm{x}}\|^{2}\}$. By $(\mu+\lambda)$-strong convexity of $g(\cdot)+\frac{\lambda}{2}\|\cdot-\hat{\bm{x}}\|^{2}$

$$ $\displaystyle g(\bm{x})+\frac{\lambda}{2}\|\bm{x}-\hat{\bm{x}}\|^{2}\geq g(\bm {x}^{\star})+\frac{\lambda}{2}\|\bm{x}^{\star}-\hat{\bm{x}}\|^{2}+\frac{\mu+ \lambda}{2}\|\bm{x}-\bm{x}^{\star}\|^{2},$ $\displaystyle g(\bm{x}^{\star})-g(\bm{x})\leq\frac{\lambda}{2}\left[\|\bm{x}- \hat{\bm{x}}\|^{2}-\|\bm{x}^{\star}-\hat{\bm{x}}\|^{2}-\|\bm{x}^{\star}-\bm{x} \|^{2}\right]-\frac{\mu}{2}\|\bm{x}-\bm{x}^{\star}\|^{2}.$ (7) $$

By the definition of $\bar{\bm{x}}$

$$ $\displaystyle\mathbb{E}\{g(\bar{\bm{x}})+\frac{\lambda}{2}\|\bar{\bm{x}}-\hat{ \bm{x}}\|^{2}\}\leq g({\bm{x}}^{\star})+\frac{\lambda}{2}\|{\bm{x}}^{\star}- \hat{\bm{x}}\|^{2}+\epsilon$ (8) $$

Combining (7) and (8) gives

$$ $\displaystyle\mathbb{E}\left\{g(\bar{\bm{x}})-g({\bm{x}})\right\}\leq\frac{ \lambda}{2}\left[\|\bm{x}-\hat{\bm{x}}\|^{2}-\|\bm{x}^{\star}-{\bm{x}}\|^{2}- \mathbb{E}\|\hat{\bm{x}}-\bar{\bm{x}}\|^{2}\right]-\frac{\mu}{2}\|\bm{x}-\bm{x }^{\star}\|^{2}+\epsilon$ (9) $\displaystyle=\mathbb{E}\left\{\frac{\lambda}{2}\left[\|\bm{x}-\hat{\bm{x}}\|^ {2}-\|\bm{x}-\bar{\bm{x}}\|^{2}-\|\bar{\bm{x}}-\hat{\bm{x}}\|^{2}\right]-\frac {\mu}{2}\|\bm{x}-\bar{\bm{x}}\|^{2}\right.$ $\displaystyle\quad\left.+\frac{\lambda+\mu}{2}\left[\|\bm{x}-\bar{\bm{x}}\|^{2 }-\|\bm{x}-\bm{x}^{\star}\|^{2}\right]+\epsilon\right\}$ $\displaystyle\leq\mathbb{E}\left\{\frac{\lambda}{2}\left[\|\bm{x}-\hat{\bm{x}} \|^{2}-\|\bm{x}-\bar{\bm{x}}\|^{2}-\|\bar{\bm{x}}-\hat{\bm{x}}\|^{2}\right]- \frac{\mu}{2}\|\bm{x}-\bar{\bm{x}}\|^{2}\right.$ $\displaystyle\quad\left.+(\lambda+\mu)\|\bm{x}-\bar{\bm{x}}\|\|\bar{\bm{x}}- \bm{x}^{\star}\|+\epsilon\right\},$ (10) $$

where the last inequality is due to the fact that $\|a\|^{2}-\|b\|^{2}\leq-2\langle a,b-a\rangle\leq 2\|a\|\|b-a\|$.

Let $\bm{x}=\bar{\bm{x}}$ in (9), and take the expectation with respect to $\bar{\bm{x}}$. Then we have

$$ $\frac{\lambda+\mu}{2}\mathbb{E}\,\|{\bm{x}}^{\star}-\bar{\bm{x}}\|^{2}\leq\epsilon.$ $$

By Hölder’s inequality we have

$$ $\displaystyle\mathbb{E}\|\bm{x}-\bar{\bm{x}}\|\|\bar{\bm{x}}-\bm{x}^{\star}\| \leq\left(\mathbb{E}\|\bm{x}-\bar{\bm{x}}\|^{2}\right)^{\frac{1}{2}}\left( \mathbb{E}\|\bm{x}^{\star}-\bar{\bm{x}}\|^{2}\right)^{\frac{1}{2}}$ $\displaystyle\leq\frac{1}{2}\left(\mathbb{E}\|\bm{x}^{\star}-\bar{\bm{x}}\|^{2 }\right)^{\frac{1}{2}}\left(1+\mathbb{E}\|\bm{x}-\bar{\bm{x}}\|^{2}\right).$ $\displaystyle\leq\frac{1}{2}\sqrt{\frac{2\epsilon}{\lambda+\mu}}\left(1+ \mathbb{E}\|\bm{x}-\bar{\bm{x}}\|^{2}\right).$ $$

Combining the above results and (10) we get the desired result.
∎

Consider solving the problem from Line 7 to Line 15 in Algorithm 1 while updating $\bm{w}$,

$$ $\displaystyle\min_{\bm{w}}P_{k}(\bm{w}):=g(\bm{w})+\bm{\lambda}_{k+1}^{\top} \bm{\ell}(\bm{w})+\frac{1}{2\eta_{k}}\|\bm{w}-\bm{w}_{k}\|^{2}$ $$

The following lemma provides the accuracy between $\bm{w}_{k+1}$ and $\arg\min_{\bm{w}}P_{k}(\bm{w})$.

###### Lemma 3 .

Let $P_{k}(\bm{w}):=\bm{\lambda}_{k+1}^{\top}\bm{\ell}(\bm{w})+g(\bm{w})+\frac{1}{2
\tau_{k}}\|\bm{w}-\bm{w}_{k}\|^{2}$. Set $\alpha<\frac{L}{4}$, $m_{k}=\Theta(\frac{L\tau_{k}}{\mu\tau_{k}+1})$ and $T_{k}=O(m_{k}\log\frac{1}{\epsilon})$ in Algorithm 1. The overall sample complexity of obtaining an $\epsilon$-approximate solution such that $\mathbb{E}P_{k}(\bm{w}^{k+1})-\min_{\bm{w}}P_{k}(\bm{w})\leq\epsilon$ is $O\left(\left(n+\frac{L\tau_{k}}{\mu\tau_{k}+1}\right)\log\frac{1}{\epsilon}\right)$. Moreover, we can set $\alpha=\frac{1}{12L}$ and $m_{k}=\frac{96L}{\mu+\tau_{k}^{-1}}+2$ in practice.

###### Proof.

First note that $\bm{\lambda}^{k+1\top}\bm{\ell}(\bm{w})$ is $L$-smooth since

$$ $\|\sum_{i=i}^{n}\lambda_{k+1,i}\ell_{i}(\bm{x})-\sum_{i=1}^{n}\lambda_{k+1,i} \ell_{i}(\bm{y})\|\leq\sum_{i=1}^{n}\lambda_{k+1,i}\|\ell_{i}(\bm{x})-\ell_{i} (\bm{y})\|\leq L\|\bm{x}-\bm{y}\|,$ $$

for $\forall\bm{x},\bm{y}\in\mathbb{R}^{d}$. In the last inequality we use $\sum_{i}\lambda_{k+1,i}\leq 1$ due to $\bm{\lambda}_{k+1}\in\Pi_{\bm{\sigma}}$ and $L$-smoothness of $\ell_{i}$. Moreover, it is not hard to see that $P_{k}(\bm{w})$ is $\mu+\tau_{k}^{-1}$-strongly convex. By [, Theorem 1] we get the desired result.
∎

Next, for convenience, we consider that $T_{k}$ (will be determined in Theorem 2) is large enough so that $\bm{w}_{k}$ is a $\delta_{k}$-near optimal solution of $P_{k}(\bm{w})$, that is, $\mathbb{E}_{k}P_{k}(\bm{w}_{k+1})-\min_{\bm{w}}P_{k}(\bm{w})\leq\epsilon$. Here, $\mathbb{E}_{k}$ represents the conditional expectation with respect to the random samples used to compute $\bm{w}_{k+1}$ given $\bm{w}_{k},\dots,\bm{w}_{0}$. Then, we can provide a one-step analysis of the outer loop of our algorithm. We use $L(\bm{w},\bm{\lambda})=\bm{\lambda}^{\top}\bm{\ell}(\bm{w})+g(\bm{w})$ in the analysis, which is natural because our algorithm iteratively updates $\bm{w}$ and $\bm{\lambda}$.

###### Lemma 4 .

Suppose Assumption 1 holds. Let $\{\bm{w}_{k}\}$ and $\{\bm{\lambda}_{k}\}$ be the sequences generated by Algorithm 1. Then for any $\bm{w}\in\mathbb{R}^{d}$ and $\bm{\lambda}\in\Pi_{\bm{\sigma}}$, the following inequality holds,

$$ $\displaystyle\mathbb{E}_{k}\,\left\{L(\bm{w}_{k+1},\bm{\lambda})-L(\bm{w},\bm{ \lambda}_{k+1})\right\}$ $\displaystyle\leq$ $\displaystyle\mathbb{E}_{k}\left\{\langle\bm{\lambda}-\bm{\lambda}_{k+1},\bm{ \ell}(\bm{w}_{k+1})\rangle+\frac{1}{2\eta_{k}}\left[\|\bm{\lambda}-\bm{\lambda }_{k}\|^{2}-\|\bm{\lambda}-\bm{\lambda}_{k+1}\|^{2}-\|\bm{\lambda}_{k+1}-\bm{ \lambda}_{k}\|^{2}\right]\right.$ $\displaystyle+\frac{1}{2\tau_{k}}\left[\|\bm{w}-\bm{w}_{k}\|^{2}-\|\bm{w}-\bm{ w}_{k+1}\|^{2}-\|\bm{w}_{k+1}-\bm{w}_{k}\|^{2}\right]-\frac{\mu}{2}\|\bm{w}- \bm{w}_{k+1}\|^{2}$ $\displaystyle+\left.\langle\bm{v}_{k},\bm{\lambda}_{k+1}-\bm{\lambda}\rangle+ \delta_{k}+\sqrt{\frac{(\tau_{k}^{-1}+\mu)\delta_{k}}{2}}(\|\bm{w}-\bm{w}_{k+1 }\|^{2}+1)\right\},$ $$

where $\mathbb{E}_{k}$ represents the conditional expectation with respect to the random samples used to compute $\bm{w}_{k+1}$ given $\bm{w}_{k},\dots,\bm{w}_{0}$.

###### Proof.

From the update of $\bm{\lambda}_{k+1}$ and Lemma 2 we have

$$ $0\leq\frac{1}{2\eta_{k}}\left[\|\bm{\lambda}-\bm{\lambda}_{k}\|^{2}-\|\bm{ \lambda}-\bm{\lambda}_{k+1}\|^{2}-\|\bm{\lambda}_{k+1}-\bm{\lambda}_{k}\|^{2} \right]+\langle\bm{v}_{k},\bm{\lambda}_{k+1}-\bm{\lambda}\rangle.$ (11) $$

From the update of $\bm{w}_{k+1}$ and Lemma 2 we have

$$ $\displaystyle\mathbb{E}_{k}\left\{g(\bm{w}_{k+1})+\langle\bm{\lambda}_{k+1}, \bm{\ell}(\bm{w}_{k+1})\rangle-g(\bm{w})-\langle\bm{\lambda}_{k+1},\bm{\ell}( \bm{w})\rangle\right\}$ $\displaystyle\leq\mathbb{E}_{k}\left\{\frac{1}{2\tau_{k}}\left[\|\bm{w}-\bm{w} _{k}\|^{2}-\|\bm{w}-\bm{w}_{k+1}\|^{2}-\|\bm{w}_{k+1}-\bm{w}_{k}\|^{2}\right]- \frac{\mu}{2}\|\bm{w}-\bm{w}_{k+1}\|^{2}\right.$ $\displaystyle\left.\quad+\delta_{k}+\sqrt{\frac{(\tau_{k}^{-1}+\mu)\delta_{k}} {2}}(\|\bm{w}-\bm{w}_{k+1}\|^{2}+1)\right\}.$ (12) $$

Taking the conditional expectation $\mathbb{E}_{k}$ of both sides of (11) and summing with (A) we obtain

$$ $\displaystyle\mathbb{E}_{k}\,\left\{L(\bm{w}_{k+1},\bm{\lambda})-L(\bm{w},\bm{ \lambda}_{k+1})\right\}$ (13) $\displaystyle=$ $\displaystyle\mathbb{E}_{k}\,\left\{g(\bm{w}_{k+1})+\langle\bm{\lambda},\bm{ \ell}(\bm{w}_{k+1})\rangle-g(\bm{w})-\langle\bm{\lambda}_{k+1},\bm{\ell}(\bm{w })\rangle\right\}$ $\displaystyle\leq$ $\displaystyle\mathbb{E}_{k}\,\left\{\langle\bm{\lambda}-\bm{\lambda}_{k+1},\bm {\ell}(\bm{w}_{k+1})\rangle+\frac{1}{2\eta_{k}}\left[\|\bm{\lambda}-\bm{ \lambda}_{k}\|^{2}-\|\bm{\lambda}-\bm{\lambda}_{k+1}\|^{2}-\|\bm{\lambda}_{k+1 }-\bm{\lambda}_{k}\|^{2}\right]\right.$ $\displaystyle+\frac{1}{2\tau_{k}}\left[\|\bm{w}-\bm{w}_{k}\|^{2}-\|\bm{w}-\bm{ w}_{k+1}\|^{2}-\|\bm{w}_{k+1}-\bm{w}_{k}\|^{2}\right]+\langle\bm{v}_{k},\bm{ \lambda}_{k+1}-\bm{\lambda}\rangle$ $\displaystyle\left.-\frac{\mu}{2}\|\bm{w}-\bm{w}_{k+1}\|^{2}+\delta_{k}+\sqrt{ \frac{(\tau_{k}^{-1}+\mu)\delta_{k}}{2}}(\|\bm{w}-\bm{w}_{k+1}\|^{2}+1)\right\}.$ $$

∎

###### Lemma 5 .

Under the same assumptions as Lemma 4, for any $\bm{w}\in\mathbb{R}^{d}$ and $\bm{\lambda}\in\Pi_{\bm{\sigma}}$, we have

$$ $\displaystyle\mathbb{E}_{k}\left\{L(\bm{w}_{k+1},\bm{\lambda})-L(\bm{w},\bm{ \lambda}_{k+1})\right\}$ (14) $\displaystyle\leq$ $\displaystyle\mathbb{E}_{k}\left\{\frac{1}{2\eta_{k}}\left[\|\bm{\lambda}-\bm{ \lambda}_{k}\|^{2}-\|\bm{\lambda}-\bm{\lambda}_{k+1}\|^{2}\right]+\frac{1}{2 \tau_{k}}\|\bm{w}-\bm{w}_{k}\|^{2}-\frac{1}{2}\left(\frac{1}{\tau_{k}}+\mu \right)\|\bm{w}-\bm{w}_{k+1}\|^{2}\right.$ $\displaystyle+\langle\bm{\ell}(\bm{w}_{k+1})-\bm{\ell}(\bm{w}_{k}),\bm{\lambda }-\bm{\lambda}_{k+1}\rangle-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})-\bm{\ell}( \bm{w}_{k-1}),\bm{\lambda}-\bm{\lambda}_{k}\rangle$ $\displaystyle-\frac{1}{2}\left[\frac{1}{\eta_{k}}-\theta_{k}\frac{G}{\alpha_{k }}\right]\|\bm{\lambda}_{k}-\bm{\lambda}_{k+1}\|^{2}-\frac{1}{2\tau_{k}}\|\bm{ w}_{k}-\bm{w}_{k+1}\|^{2}+\frac{G\theta_{k}\alpha_{k}}{2}\|\bm{w}_{k}-\bm{w}_{ k-1}\|^{2}$ $\displaystyle\left.+\delta_{k}+\sqrt{\frac{(\tau_{k}^{-1}+\mu)\delta_{k}}{2}}( \|\bm{w}-\bm{w}_{k+1}\|^{2}+1)\right\}.$ $$

###### Proof.

First, we have

$$ $\displaystyle\langle\bm{v}_{k},\bm{\lambda}_{k+1}-\bm{\lambda}\rangle$ $\displaystyle=$ $\displaystyle\langle\bm{\ell}(\bm{w}_{k})+\theta_{k}\left(\bm{\ell}(\bm{w}_{k} )-\bm{\ell}(\bm{w}_{k-1})\right),\bm{\lambda}_{k+1}-\bm{\lambda}\rangle$ $\displaystyle=$ $\displaystyle-\langle\bm{\ell}(\bm{w}_{k+1}),\bm{\lambda}-\bm{\lambda}_{k+1} \rangle+\langle\bm{\ell}(\bm{w}_{k+1})-\bm{\ell}(\bm{w}_{k}),\bm{\lambda}-\bm{ \lambda}_{k+1}\rangle$ $\displaystyle-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})-\bm{\ell}(\bm{w}_{k-1}), \bm{\lambda}-\bm{\lambda}_{k}\rangle-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})- \bm{\ell}(\bm{w}_{k-1}),\bm{\lambda}_{k}-\bm{\lambda}_{k+1}\rangle.$ $$

Then we obtain that

$$ $\displaystyle\langle\bm{\lambda}-\bm{\lambda}_{k+1},\bm{\ell}(\bm{w}_{k+1}) \rangle+\langle\bm{v}_{k},\bm{\lambda}_{k+1}-\bm{\lambda}\rangle$ (15) $\displaystyle\leq\langle\bm{\ell}(\bm{w}_{k+1})-\bm{\ell}(\bm{w}_{k}),\bm{ \lambda}-\bm{\lambda}_{k+1}\rangle-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})-\bm{ \ell}(\bm{w}_{k-1}),\bm{\lambda}-\bm{\lambda}_{k}\rangle$ $\displaystyle\quad-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})-\bm{\ell}(\bm{w}_{k- 1}),\bm{\lambda}_{k}-\bm{\lambda}_{k+1}\rangle.$ $$

Next we bound the last term on the right-hand side of (15).

$$ $\displaystyle\langle\bm{\ell}(\bm{w}_{k})-\bm{\ell}(\bm{w}_{k-1}),\bm{\lambda} _{k}-\bm{\lambda}_{k+1}\rangle$ (16) $\displaystyle\leq$ $\displaystyle G\|\bm{w}_{k}-\bm{w}_{k-1}\|\|\bm{\lambda}_{k}-\bm{\lambda}_{k+1}\|$ $\displaystyle\leq$ $\displaystyle\frac{G\alpha_{k}}{2}\|\bm{w}_{k}-\bm{w}_{k-1}\|^{2}+\frac{G}{2 \alpha_{k}}\|\bm{\lambda}_{k}-\bm{\lambda}_{k+1}\|^{2},$ $$

where the first inequality is due to the Lipschitz conditions of $\ell_{i}$ and in the second inequality we use Young’s inequality with $\alpha_{k}>0$.

Combing (15) and (16) we get

$$ $\displaystyle\langle\bm{\lambda}-\bm{\lambda}_{k+1},\bm{\ell}(\bm{w}_{k+1}) \rangle+\langle\bm{v}_{k},\bm{\lambda}_{k+1}-\bm{\lambda}\rangle$ (17) $\displaystyle\leq\langle\bm{\ell}(\bm{w}_{k+1})-\bm{\ell}(\bm{w}_{k}),\bm{ \lambda}-\bm{\lambda}_{k+1}\rangle-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})-\bm{ \ell}(\bm{w}_{k-1}),\bm{\lambda}-\bm{\lambda}_{k}\rangle$ $\displaystyle\quad+\frac{G\alpha_{k}\theta_{k}}{2}\|\bm{w}_{k}-\bm{w}_{k-1}\|^ {2}+\frac{G\theta_{k}}{2\alpha_{k}}\|\bm{\lambda}_{k}-\bm{\lambda}_{k+1}\|^{2}.$ $$

Taking the conditional expectation $\mathbb{E}_{k}$ of both sides of (17) and combing it with Lemma 4
we get the desired result.
∎

We remark that $\alpha_{k}$ does not need to be computed in the actual algorithm but only exists in the theoretical analysis. Next, we try to telescope the terms on the right hand side of (14) by multiplying each term by $\gamma_{k}$. To ensure that the adjacent terms in the sequence $k=0,\dots,K-1$ can be canceled out during summation, we need the parameters of the algorithm to satisfy the following conditions.

###### Condition 1 .

For $k=0,1,...$, the following conditions for parameters in the analysis and Algorithm 1:

$$ $\displaystyle\frac{\gamma_{k+1}}{\eta_{k+1}}\leq\frac{\gamma_{k}}{\eta_{k}},$ (18a) $\displaystyle\frac{\gamma_{k+1}}{\tau_{k+1}}\leq\gamma_{k}\left(\frac{1}{\tau_ {k}}+\mu-\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right),$ (18b) $\displaystyle\gamma_{k}=\gamma_{k+1}\theta_{k+1},$ (18c) $\displaystyle G\alpha_{k+1}\leq\frac{1}{\tau_{k}},$ (18d) $\displaystyle\theta_{k}\frac{G}{\alpha_{k}}\leq\frac{1}{\eta_{k}}.$ (18e) $$

###### Lemma 6 .

Assume Assumption 1 holds and Condition 1 is satisfied. Then for all $\bm{w}\in\mathbb{R}^{d}$ and $\bm{\lambda}\in\Pi_{\bm{\sigma}}$ we have

$$ $\frac{\gamma_{K}}{2\tau_{K}}\mathbb{E}\|\bm{w}^{\star}-\bm{w}_{K}\|^{2}\leq \frac{\gamma_{0}}{2\eta_{0}}\|\bm{\lambda}^{\star}-\bm{\lambda}_{0}\|^{2}+ \frac{\gamma_{0}}{2\tau_{0}}\|\bm{w}^{\star}-\bm{w}_{0}\|^{2}+\sum_{k=0}^{K-1} \left(\delta_{k}\gamma_{k}+\frac{\gamma_{k}}{2}\sqrt{2(\mu+\tau_{k}^{-1}) \delta_{k}}\right),$ $$

where $\bm{w}^{\star}=\arg\min_{\bm{w}}R_{\bm{\sigma}}(\bm{w})+g(\bm{w})$ and $\bm{\lambda}^{\star}=\bm{\sigma}_{\pi^{-1}}$. $\pi$ is the permutation that arranges $\ell_{1}(\bm{w}^{\star}),\dots,\ell_{n}(\bm{w}^{\star})$ in ascending order, that is, $\ell_{\pi(1)}(\bm{w}^{\star})\leq\dots\leq\ell_{\pi(n)}(\bm{w}^{\star})$.

###### Proof.

Taking expectations with respect to $\bm{w}_{k},\dots,\bm{w}_{1}$
in (14) and using the law
of total expectation yields

$$ $\displaystyle\mathbb{E}\left\{L(\bm{w}_{k+1},\bm{\lambda})-L(\bm{w},\bm{ \lambda}_{k+1})\right\}$ (19) $\displaystyle\leq$ $\displaystyle\mathbb{E}\,\left\{\frac{1}{2\eta_{k}}\left[\|\bm{\lambda}-\bm{ \lambda}_{k}\|^{2}-\|\bm{\lambda}-\bm{\lambda}_{k+1}\|^{2}\right]-\frac{1}{2} \left[\frac{1}{\eta_{k}}-\frac{G\theta_{k}}{\alpha_{k}}\right]\|\bm{\lambda}_{ k}-\bm{\lambda}_{k+1}\|^{2}\right.$ $\displaystyle+\langle\bm{\ell}(\bm{w}_{k+1})-\bm{\ell}(\bm{w}_{k}),\bm{\lambda }-\bm{\lambda}_{k+1}\rangle-\theta_{k}\langle\bm{\ell}(\bm{w}_{k})-\bm{\ell}( \bm{w}_{k-1}),\bm{\lambda}-\bm{\lambda}_{k}\rangle$ $\displaystyle+\frac{1}{2\tau_{k}}\|\bm{w}-\bm{w}_{k}\|^{2}-\frac{1}{2}\left( \frac{1}{\tau_{k}}+\mu-\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right)\|\bm{w}- \bm{w}_{k+1}\|^{2}$ $\displaystyle\left.-\frac{1}{2\tau_{k}}\|\bm{w}_{k}-\bm{w}_{k+1}\|^{2}+\frac{G \theta_{k}\alpha_{k}}{2}\|\bm{w}_{k}-\bm{w}_{k-1}\|^{2}+\delta_{k}+\frac{1}{2} \sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right\}\,.$ $$

Multiplying both sides of (19) by $\gamma_{k}$ and summing over $k=0$ to $K-1$ we obtain that

$$ $\displaystyle\sum_{k=0}^{K-1}\gamma_{k}\mathbb{E}\left\{L(\bm{w}_{k+1},\bm{ \lambda})-L(\bm{w},\bm{\lambda}_{k+1})\right\}$ $\displaystyle\leq$ $\displaystyle\mathbb{E}\left\{\ \frac{\gamma_{0}}{2\eta_{0}}\|\bm{\lambda}-\bm {\lambda}_{0}\|^{2}+\sum_{k=0}^{K-2}\underbrace{\frac{1}{2}\left(\frac{\gamma_ {k+1}}{\eta_{k+1}}-\frac{\gamma_{k}}{\eta_{k}}\right)\|\bm{\lambda}-\bm{ \lambda}_{k+1}\|^{2}}_{A}-\frac{\gamma_{K-1}}{2\eta_{K-1}}\|\bm{\lambda}-\bm{ \lambda}_{K}\|^{2}\right.$ $\displaystyle+$ $\displaystyle\frac{\gamma_{0}}{2\tau_{0}}\|\bm{w}-\bm{w}_{0}\|^{2}+\sum_{k=0}^ {K-2}\underbrace{\frac{1}{2}\left[\frac{\gamma_{k+1}}{\tau_{k+1}}-\gamma_{k} \left(\frac{1}{\tau_{k}}+\mu-\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right) \right]\|\bm{w}-\bm{w}_{k+1}\|^{2}}_{B}$ $\displaystyle-$ $\displaystyle\frac{\gamma_{K-1}}{2}\left(\frac{1}{\tau_{K-1}}+\mu-\sqrt{2(\mu+ \tau_{K-1}^{-1})\delta_{K-1}}\right)\|\bm{w}-\bm{w}_{K}\|^{2}$ $\displaystyle+$ $\displaystyle\sum_{k=0}^{K-2}\underbrace{(\gamma_{k}-\gamma_{k+1}\theta_{k+1}) }_{C}\langle\bm{\ell}(\bm{w}_{k+1})-\bm{\ell}(\bm{w}_{k}),\bm{\lambda}-\bm{ \lambda}_{k+1}\rangle+\gamma_{K-1}\langle\bm{\ell}(\bm{w}_{K})-\bm{\ell}(\bm{w }_{K-1}),\bm{\lambda}-\bm{\lambda}_{K}\rangle$ $\displaystyle+$ $\displaystyle\frac{1}{2}\sum_{k=0}^{K-2}\underbrace{\left(\gamma_{k+1}\theta_{ k+1}\alpha_{k+1}G-\frac{\gamma_{k}}{\tau_{k}}\right)}_{D}\|\bm{w}_{k}-\bm{w}_{ k+1}\|^{2}-\frac{\gamma_{K-1}}{2\tau_{K-1}}\|\bm{w}_{K}-\bm{w}_{K-1}\|^{2}$ $\displaystyle+$ $\displaystyle\left.\frac{1}{2}\sum_{k=0}^{K-1}\underbrace{\left[-\gamma_{k} \left(\frac{1}{\eta_{k}}-\theta_{k}\frac{G}{\alpha_{k}}\right)\right]}_{E}\| \bm{\lambda}_{k}-\bm{\lambda}_{k+1}\|^{2}+\sum_{k=0}^{K-1}\left(\delta_{k} \gamma_{k}+\frac{\gamma_{k}}{2}\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right) \right\}.$ $$

Here we use $\bm{\ell}(\bm{w}_{0})-\bm{\ell}(\bm{w}_{-1})=0$ by $\bm{w}_{0}=\bm{w}_{-1}$ and $\bm{\lambda}_{0}=\bm{\lambda}_{-1}$. By Condition 1, we have $A,B,D,E\leq 0$ and $C=0$.

Then we have

$$ $\displaystyle\sum_{k=0}^{K-1}\gamma_{k}\mathbb{E}\left\{L(\bm{w}_{k+1},\bm{ \lambda})-L(\bm{w},\bm{\lambda}_{k+1})\right\}$ (20) $\displaystyle\leq$ $\displaystyle\mathbb{E}\left\{\frac{\gamma_{0}}{2\eta_{0}}\|\bm{\lambda}-\bm{ \lambda}_{0}\|^{2}-\frac{\gamma_{K-1}}{2\eta_{K-1}}\|\bm{\lambda}-\bm{\lambda} _{K}\|^{2}+\frac{\gamma_{0}}{2\tau_{0}}\|\bm{w}-\bm{w}_{0}\|^{2}\right.$ $\displaystyle-$ $\displaystyle\frac{\gamma_{K-1}}{2}\left(\frac{1}{\tau_{K-1}}+\mu-\sqrt{2(\mu+ \tau_{K-1}^{-1})\delta_{K-1}}\right)\|\bm{w}-\bm{w}_{K}\|^{2}$ $\displaystyle+$ $\displaystyle\gamma_{K-1}\langle\bm{\ell}(\bm{w}_{K})-\bm{\ell}(\bm{w}_{K-1}), \bm{\lambda}-\bm{\lambda}_{K}\rangle-\frac{\gamma_{K-1}}{2\tau_{K-1}}\|\bm{w}_ {K}-\bm{w}_{K-1}\|^{2}$ $\displaystyle+$ $\displaystyle\left.\sum_{k=0}^{K-1}\left(\delta_{k}\gamma_{k}+\frac{\gamma_{k} }{2}\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right)\right\}.$ $$

Next we bound $\gamma_{K-1}\langle\bm{\ell}(\bm{w}_{K})-\bm{\ell}(\bm{w}_{K-1}),\bm{\lambda}-
\bm{\lambda}_{K}\rangle$ similar to (16). We have

$$ $\displaystyle\langle\bm{\ell}(\bm{w}_{K})-\bm{\ell}(\bm{w}_{K-1}),\bm{\lambda} -\bm{\lambda}_{K}\rangle\leq\frac{G\alpha_{K}}{2}\|\bm{w}_{K}-\bm{w}_{K-1}\|^{ 2}+\frac{1}{2}\frac{G}{\alpha_{K}}\|\bm{\lambda}-\bm{\lambda}_{K}\|^{2}.$ $$

Taking the expectation and plugging this into (20) we obtain that

$$ $\displaystyle\sum_{k=0}^{K-1}\gamma_{k}\mathbb{E}\,\left\{L(\bm{w}_{k+1},\bm{ \lambda})-L(\bm{w},\bm{\lambda}_{k+1})\right\}$ (21) $\displaystyle\leq$ $\displaystyle\mathbb{E}\,\left\{\frac{\gamma_{0}}{2\eta_{0}}\|\bm{\lambda}-\bm {\lambda}_{0}\|^{2}-\frac{1}{2}\underbrace{\left[\frac{\gamma_{K-1}}{\eta_{K-1 }}-\gamma_{K-1}\frac{G}{\alpha_{K}}\right]}_{\tilde{A}}\|\bm{\lambda}-\bm{ \lambda}_{K}\|^{2}\right.$ $\displaystyle+$ $\displaystyle\frac{\gamma_{0}}{2\tau_{0}}\|\bm{w}-\bm{w}_{0}\|^{2}-\underbrace {\frac{\gamma_{K-1}}{2}\left(\frac{1}{\tau_{K-1}}+\mu-\sqrt{2(\mu+\tau_{K-1}^{ -1})\delta_{K-1}}\right)}_{\tilde{B}}\|\bm{w}-\bm{w}_{K}\|^{2}$ $\displaystyle+$ $\displaystyle\left.\frac{\gamma_{K-1}}{2}\underbrace{\left(\alpha_{K}G-\frac{1 }{\tau_{K-1}}\right)}_{\tilde{C}}\|\bm{w}_{K}-\bm{w}_{K-1}\|^{2}+\sum_{k=0}^{K -1}\left(\delta_{k}\gamma_{k}+\frac{\gamma_{k}}{2}\sqrt{2(\mu+\tau_{k}^{-1}) \delta_{k}}\right)\right\}$ $$

We analyze $\tilde{A}$-$\tilde{D}$ under Condition 1:

$$ $\displaystyle\tilde{A}\overset{(\ref{assumption:step-size_1})}{\geq}\left[ \frac{\gamma_{K}}{\eta_{K}}-\gamma_{K-1}\frac{G}{\alpha_{K}}\right]\overset{( \ref{assumption:step-size_3})}{=}\gamma_{K}\left[\frac{1}{\eta_{K}}-\theta_{K} \frac{G}{\alpha_{K}}\right]\overset{(\ref{assumption:step-size_7})}{\geq}0,$ $\displaystyle\tilde{B}\overset{(\ref{assumption:step-size_2})}{\geq}\frac{ \gamma_{K}}{2\tau_{K}},$ $\displaystyle\tilde{C}\overset{(\ref{assumption:step-size_6})}{\leq}0.$ $$

We obtain that

$$ $\displaystyle\sum_{k=0}^{K-1}\gamma_{k}\mathbb{E}\,\left\{L(\bm{w}_{k+1},\bm{ \lambda})-L(\bm{w},\bm{\lambda}_{k+1})\right\}$ $\displaystyle\leq$ $\displaystyle\frac{\gamma_{0}}{2\eta_{0}}\|\bm{\lambda}-\bm{\lambda}_{0}\|^{2} +\frac{\gamma_{0}}{2\tau_{0}}\|\bm{w}-\bm{w}_{0}\|^{2}-\frac{\gamma_{K}}{2\tau _{K}}\mathbb{E}\,\|\bm{w}-\bm{w}_{K}\|^{2}+\sum_{k=0}^{K-1}\left(\delta_{k} \gamma_{k}+\frac{\gamma_{k}}{2}\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\right)$ $$

Let $\bm{w}=\bm{w}^{\star},\bm{\lambda}=\bm{\lambda}^{\star}$, for any $\bm{w}\in\mathbb{R}^{d}$ and $\bm{\lambda}\in\Pi_{\bm{\sigma}}$, we have $L(\bm{w}^{\star},\bm{\lambda}^{\star})=\max_{\bm{\lambda}\in\Pi_{\bm{\sigma}}}
L(\bm{w}^{\star},\bm{\lambda})\geq L(\bm{w}^{\star},\bm{\lambda})$. On the other hand, we have $L(\bm{w},\bm{\lambda}^{\star})\geq L(\bm{w}^{\star},\bm{\lambda}^{\star})=\min
_{\bm{w}}L(\bm{w},\bm{\lambda}^{\star})$. Thus we obtain that $L(\bm{w}_{k+1},\bm{\lambda})-L(\bm{w},\bm{\lambda}_{k+1})\geq 0$ for $\forall k=0,\dots,K-1$. This completes the proof.
∎

Now we are ready to present our main theorem.
By choosing appropriate parameters in Algorithm 1 to satisfy Condition 1, we can achieve the desired convergence rate.

###### Theorem 2 .

Assume Assumption 1 holds.
Set $\eta_{k}=\frac{\mu(k+1)}{8G^{2}},\theta_{k}=\frac{k}{k+1},\gamma_{k}=k+1,\tau_
{k}=\frac{4}{\mu(k+1)},\alpha_{k+1}=G\eta_{k},\delta_{k}=\min\left(\frac{\mu}{
8(k+5)},\mu(k+1)^{-6}\right)$.
Then Condition 1 is satisfied by the above parameters. Moreover, use the same conditions as Lemma 3, while replacing $\epsilon$ with $\delta_{k}$ for $k=0,1,\dots,K-1$. Then we have

$$ $\mathbb{E}\|\bm{w}_{K}-\bm{w}^{\star}\|^{2}=O\left(\frac{G^{2}}{\mu^{2}K^{2}} \right).$ $$

###### Proof.

First, we obtain an $\delta_{k}$ approximate solution to (A) through $T_{k}$ updates to $\bm{w}$ in Algorithm 1 by Lemma 3.
We then verify that Condition 1 is satisfied by the parameters.

It is not hard to see that $\frac{\gamma_{k+1}}{\gamma_{k}}=\frac{\eta_{k+1}}{\eta_{k}}=\frac{k+2}{k+1}$ and $\theta_{k+1}=\frac{\gamma_{k}}{\gamma_{k+1}}=\frac{k+1}{k+2}$. Thus (18a) and (18c) are satisfied.

Since $\delta_{k}\leq\frac{\mu}{8(k+5)}$, we have $\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}=\sqrt{2\mu(1+\frac{k+1}{4})\delta_{k}}
\leq\frac{\mu}{4}$. Then we obtain that

$$ $\frac{\gamma_{k+1}}{\gamma_{k}\tau_{k+1}}=\frac{k+2}{4}\mu+\frac{k+2}{4(k+1)} \mu\leq\frac{k+4}{4}\mu,$ $$

and

$$ $\frac{1}{\tau_{k}}+\mu-\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}}\geq\frac{k+1}{4} \mu+\mu-\frac{\mu}{4}=\frac{k+4}{4}\mu.$ $$

Thus (18b) holds.

Furthermore, (18d) and (18e) hold due to $G\alpha_{k+1}=G^{2}\eta_{k}=\frac{k+1}{8}\mu\leq\frac{k+1}{4}\mu=\frac{1}{\tau
_{k}}$ and $\theta_{k}\frac{G}{\alpha_{k}}=\frac{\eta_{k-1}}{\eta_{k}}\frac{G}{G\eta_{k-1}
}=\frac{1}{\eta_{k}}$.

Now Condition 1 is satisfied. By Lemma 6, we have

$$ $\frac{\gamma_{K}}{2\tau_{K}}\mathbb{E}\,\|\bm{w}-\bm{w}_{K}\|^{2}\leq\frac{ \gamma_{0}}{2\eta_{0}}\|\bm{\lambda}^{\star}-\bm{\lambda}_{0}\|^{2}+\frac{ \gamma_{0}}{2\tau_{0}}\|\bm{w}^{\star}-\bm{w}_{0}\|^{2}+\sum_{k=0}^{K-1}\left( \delta_{k}\gamma_{k}+\frac{\gamma_{k}}{2}\sqrt{2(\mu+\tau_{k}^{-1})\delta_{k}} \right).$ $$

Since $\delta_{k}\leq\mu(k+1)^{-6}$, we have $\sum_{k=0}^{\infty}\delta_{k}\gamma_{k}\leq\mu\sum_{k=0}^{\infty}(k+1)^{-5}
\leq\frac{\mu}{4}$,
and

$$ $\sum_{k=0}^{\infty}\gamma_{k}\sqrt{(\mu+\tau_{k}^{-1})\delta_{k}}\leq\frac{ \sqrt{2}\mu}{4}\sum_{k=0}^{\infty}(k+1)^{-2}\sqrt{k+5}\leq\frac{\sqrt{2}\mu}{4 }\sum_{k=0}^{\infty}(k+1)^{-2}\left(\sqrt{k+1}+2\right)\leq\sqrt{2}\mu.$ $$

Finally, by $\frac{\gamma_{K}}{2\tau_{K}}=\frac{\mu(K+1)^{2}}{8},\tau_{0}=\frac{4}{\mu}$ and $\eta_{0}=\frac{\mu}{8G^{2}}$, we get the desired result.
∎

###### Corollary 2 .

Suppose the step-size $\alpha<\frac{L}{4}$ and $m_{k}=\Theta(\frac{L}{\mu(k+1)})$. With the same conditions as Theorem 2, we obtain an output $\bm{w}_{K}$ such that $\mathbb{E}\|\bm{w}_{K}-\bm{w}^{\star}\|^{2}\leq\epsilon$ in a total sample complexity of $O\left(n\frac{G}{\mu\sqrt{\epsilon}}\log\frac{G}{\mu^{2}\sqrt{\epsilon}}+\frac
{L}{\mu}\log\frac{G}{\mu\sqrt{\epsilon}}\log\frac{G}{\mu^{2}\sqrt{\epsilon}}\right)$ using Algorithm 1.

###### Proof.

Recall that $\tau_{k}=\frac{4}{\mu(k+1)}$. It is not hard to see that $\frac{L\tau_{k}}{\mu\tau_{k}+1}=\frac{4L}{\mu(k+5)}\leq\frac{4L}{\mu(k+1)}$. By Lemma 3, we get a $\delta_{k}$ approximate solution with the sample complexity of $C_{\bm{w}_{k+1}}=O\left(\left(n+\frac{L}{\mu(k+1)}\right)\log\left(\delta_{k}^
{-1}\right)\right)$. We set $\delta_{k}=\min\left(\frac{\mu}{8(k+5)},\mu(k+1)^{-6}\right)=\mu(k+1)^{-6}$ for $k\geq 1$. And $\delta_{0}=\mu/40$. The total sample complexity is

$$ $\displaystyle\sum_{k=0}^{K-1}C_{\bm{w}_{k+1}}=\sum_{k=0}^{K-1}O\left(\left(n+ \frac{L}{\mu(k+1)}\right)\left(\log(k+1)-\log\mu\right)\right)$ $\displaystyle=O\left(nK\log\frac{K}{\mu}+\frac{L}{\mu}\log K\log\frac{K}{\mu} \right).$ $$

In the last equality, we calculate $\sum_{k=1}^{K}\frac{\log k}{k}=O\left(\left(\log K\right)^{2}\right)$, $\sum_{k=1}^{K}\log k=O\left(K\log K\right)$ and $\sum_{k=1}^{K}\frac{1}{k}=O(\log K)$. By Theorem 2, to achieve an $\epsilon$ approximate solution, we need $K=O\left(\frac{G}{\mu\sqrt{\epsilon}}\right)$. Therefore, the total sample complexity is $O\left(n\frac{G}{\mu\sqrt{\epsilon}}\log\frac{G}{\mu^{2}\sqrt{\epsilon}}+\frac
{L}{\mu}\log\frac{G}{\mu\sqrt{\epsilon}}\log\frac{G}{\mu^{2}\sqrt{\epsilon}}
\right).$
∎

###### Lemma 2 .

###### Proof.

###### Lemma 3 .

###### Proof.

###### Lemma 4 .

###### Proof.

###### Lemma 5 .

###### Proof.

###### Condition 1 .

###### Lemma 6 .

###### Proof.

###### Theorem 2 .

###### Proof.

###### Corollary 2 .

###### Proof.

## Appendix B Experimental Details

We now outline the details of our experimental setup. Our experimental setup mainly follows that of Mehta et al. [].

#### Datasets.

We use the same five datasets from the regression task in Mehta et al. []. The statistical characteristics are summarized in Table 2. Tasks for each dataset are as follows:

- 1.
yacht: predicting the residual resistance of sailing yachts based on their physical attributes.
- 2.
energy: predicting the cooling load of buildings based on their physical attributes.
- 3.
concrete: predicting the compressive strength of concrete types based on their physical and chemical properties.
- 4.
kin8nm: predicting the distance of a robot arm, consisting of 8 fully rotating joints, to a target point in space.
- 5.
power: predicting the net hourly electrical energy output of a power plant based on environmental factors.

**Table 2: Statistical details of five real datasets and sources.**
| Dataset | # features | # samples | Source |
| --- | --- | --- | --- |
| yacht | 6 | 244 | Tsanas and Xifara [] |
| energy | 8 | 614 | Baressi Šegota et al. [] |
| concrete | 8 | 824 | Yeh [] |
| kin8nm | 8 | 6,553 | Akujuobi and Zhang [] |
| power | 4 | 7,654 | Tüfekci [] |

In the experiments, the sample matrix $\bm{X}\in\mathbb{R}^{n\times d}$ of the training set is standardized such that each column has a zero mean and a unit variance.

#### Objectives.

We use linear models in our experiments. For spectral risks, we adopt three types: ESRM ($\rho=2$), Extremile ($r=2.5$), and CVaR ($\alpha=0.5$), as specified in Table 1. Additionally, we set the regularizer $g(\bm{w})$ to $\frac{\mu}{2}\|\bm{w}\|^{2}$ with $\mu=\frac{1}{n}$. Thus, Problem 1 can be written as

$$ $\min_{\bm{w}}\sum_{i=1}^{n}\sigma_{i}\ell_{[i]}(\bm{w})+\frac{\mu}{2}\|\bm{w} \|^{2},$ $$

where $\ell_{i}(\bm{w})=\frac{1}{2}\|y_{i}-\bm{w}^{\top}\bm{x}_{i}\|^{2}$ and $(\bm{x}_{i},y_{i})$ are samples from the training set.

#### Hyperparameter Selection.

We use the same hyperparameter selection method as in []. We set the batch size for SGD to $64$. Each time, we randomly select a minibatch indexes from $\{1,\dots,n\}$ without replacement for SGD and randomly select a single index for other methods. For the selection of step size $\alpha$, we set the random seed $s\in\{1,\dots,S\}$. For a single seed $s$, we calculate the average training loss of the last ten epochs, donated by $L_{s}(\alpha)$. We choose $\alpha$ that minimizes $\frac{1}{S}\sum_{s=1}^{S}L_{s}(\alpha)$, where $\alpha\in\{1\times 10^{-4},3\times 10^{-4},1\times 10^{-3},3\times 10^{-3},1
\times 10^{-2},3\times 10^{-2},1\times 10^{-1},3\times 10^{-1}\}$. For LSVRG, we set the length of an epoch to $n$. For SOREL,
we set $T_{k}=m_{k}=n$. Moreover, we set batch size to 64 for all algorithms with minibatching.

For SOREL, we follow the parameter values given in Theorem 1. In particular, we set $\theta_{k}=\frac{k}{k+1}$ and $\tau_{k}=\frac{20n}{k+1}$ in all experiments. Therefore, there are only two parameters $\alpha$ and $\eta_{k}$ left to tune. We set $\eta_{k}=\frac{C(k+1)}{n}$ and choose $C$ from $\{1\times 10^{-2},2\times 10^{-2},4\times 10^{-2},1\times 10^{-1},2\times 10^{
-1},4\times 10^{-1},1\times 10^{0},2\times 10^{0},4\times 10^{0}\}$, since the Lipschitz constant $G$ is hard to estimate. We use grid search to select $\alpha$ and $C$, with the selection criteria being the same as the previous paragraph.

#### Experimental Environment.

We run all experiments on a laptop with 16.0 GB RAM and Intel i7-1360P 2.20 GHz CPU. All algorithms are implemented in Python 3.8.

## Appendix C Examples

To illustrate the necessity of stabilizing the trajectory of the primal variable in Section 3.1, we provide a toy example. For simplicity, we consider a one-dimensional problem

$$ $\min_{w\in\mathbb{R}}\sigma_{1}\ell_{[1]}(w)+\sigma_{2}\ell_{[2]}(w),$ (22) $$

where $\sigma_{1}=0,\sigma_{2}=1$ and $\ell_{1}=\frac{1}{2}(w-1)^{2},\ell_{2}(w)=\frac{1}{2}(w+1)^{2}$.
We use the following deterministic method, similar to Algorithm 1.

###### Example 1 .

For any $0<\alpha<2$, suppose we solve Problem (22) using Algorithm 2
and $T$ is sufficiently large. In that case, the iterative sequence $\{w_{k}\}$ can not converge to the optimal solution for any initial point $w_{0}$.

Figure: Algorithm 2 Simplified Algorithm for Solving the Example Problem.

Without loss of generality, we assume $w_{0}>0$, in which case $\bm{\lambda}_{1}=[0,1]^{\top}$. We solve $\min_{w_{1}\in\mathbb{R}}\frac{1}{2}(w_{1}+1)^{2}$ through sufficient steps of gradient descent to obtain $w_{1}=-1$. At this point, $\bm{\lambda}_{2}=[1,0]^{\top}$. By iterating this process, $w_{k}$ always oscillates between -1 and 1, unable to converge to $w^{\star}=0$. If $w_{0}=0$, we set $\sigma_{1}=1$ and $\sigma_{2}=0$, reaching the same conclusion. A similar conclusion can be extended to stochastic methods in the expectation sense.

We know that $\ell_{1}(w^{\star})=\ell_{2}(w^{\star})$ at the optimal point $w^{\star}=0$. Clearly, the iterative sequence of the algorithm oscillates at $w^{\star}$ and cannot converge to the optimal solution.
Although Mehta et al. [, ] employ a similar approach to update $\bm{\lambda}$ for subgradient estimations, they consider the smoothed spectral risk by adding a strongly concave term with respect to $\bm{\lambda}$. However, for the original spectral risk minimization problems, updating $\bm{\lambda}$ with their method results in discontinuities, thereby lacking convergence guarantees.

###### Example 1 .

## Peer Reviews

confidence: 3
contribution: 2
initial_rating: 6
presentation: 3
questions: See the section above.
soundness: 3
strengths: This paper is well-written and easy to follow. The claimed contribution is interesting to the community. However, I have a major concern regarding its correctness which I elaborate on in the next section.
summary: This paper proposes a stochastic algorithm for spectral risk minimization with trajectory stabilization for the primal variable. It is claimed that their approach enjoys a near-optimal rate of convergence that matches the known lower bound. They also demonstrate that this method outperforms existing baselines in various experimental settings.
weaknesses: I didn't read the proofs closely but the main theorem seems to be incorrect by a sanity check. To be more specific, $\mu$ is in unit loss/par^2 and $\delta_k$ is in unit loss, but $\delta_k \sim \mu$ in Theorem 1. Moreover, $G$ is in unit loss/par and $\epsilon$ is in unit par^2, so $G/(\mu \sqrt{\epsilon})$ is unitless. However, there is $\log{(G/(\mu^2 \sqrt{\epsilon}))}$ in the sample complexity given in Cor. 1 which makes it not unitless.

Other comments:
1. I do not see error bars (standard errors) reported in all the experiments. Are these results based on a single run or multiple independent runs? Can you report error bars so that we can tell how significant the improvements are?
2. $w^*$ is not defined in Theorem 1.

---

confidence: 3
contribution: 3
initial_rating: 6
presentation: 3
questions: 1. Can you define `total sample complexity' in Corollary 1 please?  Is this the sum of calls to a stochastic gradient oracle?  It would be helpful if this were made rigorous as my understanding of the key contribution of the present work is that the dependence in complexity on the number of samples is reduced substantially.

2. See weaknesses 1 and 2.
soundness: 3
strengths: This paper sets out an interesting problem that has some (at least distant) applications to relevant areas in ML including risk sensitivity, fairness, and generalization OOD.  The proposed solution of reducing to a minimax game and then running a primal-dual algorithm is interesting and the empirical validation uses a reasonable diversity of datasets and loss functions that suffices to convince me that the proposed algorithm is superior to alternatives, at least when considering ridge regression.
summary: This paper considers the problem of optimization in *spectral risk minimization*, which aims to produce more risk-sensitive models than standard ERM.  More precisely, spectral risk is defined as a weighted empirical risk minimization problem, where the weights depend on the order statistics of a given empirical sample, thereby subsuming both empirical risk (with all weights being the same) and such risk-sensitive measures as CVaR and Extremiles.  While minimizing spectral risk can be attractive due to the risk-sensitivity, it can be challenging to scale as the evaluation of the spectral risk of a given parameter depends on the entire dataset through the order statistics and is thus not amenable to naive applications of existing SGD approaches.  In particular, in order to obtain an unbiased estimate of the subgradient of the loss, all samples must have losses evaluated.  The authors circumvent this problem by applying the rearrangement inequality to reduce the problem to finding the value of a min-max game and then propose applying what amounts to a regularized form of projected gradient ascent-descent with momentum.  The authors present rates on their algorithm that are close to earlier lower bounds and conclude the paper with a number of experiments with ridge regression on a variety of instantiations of spectral risks and datasets, demonstrating that their approach is generally superior to alternatives.
weaknesses: There are several main weaknesses with the paper:

1. With respect to the theory, I think more discussion of the strongly convex regularizer is necessary.  I understand that it is necessary in order to ensure identifiability, but I find the discussion under Corollary 1 confusing.  For example, in the comparison to related work, the authors cite several works that set the regularization term to be $O(\epsilon)$ in order to be $\epsilon$-suboptimal with respect to the unregularized solution and note that this leads to a worse sample complexity of $O(1/\epsilon)$; I am confused as to why the current proposal does not suffer the same drawback.  Indeed, in order to be $\epsilon$-suboptimal with respect to the unperturbed spectral risk, it seeems as if the bound in Corollary 1 would scale with $\epsilon^{-3/2}$ which is obviously worse than the alternative.

2. I think a more thorough discussion of the computational complexity of SOREL is in order.  While I understand that for most steps, the gradient updates are cheap and independent of $n$, lines 2, 7, and 11 of Algorithm 1 require linear in the number of samples time.  The authors claim that these steps are for variance reduction, but do not investigate the necessity thereof.  Another way of saying this is that the authors claim almost optimality of their algorithm in how the iteration complexity scales with $\epsilon$ but it seems like the dependence on $n$ is also fairly important (and is what motivates the study itself) and the extent to which the authors' algorithm is optimal in this respect is left largely undiscussed.

3. The experiments are all restricted to linear models.  While I understand that the authors wish to ensure that their key Assumption 1 is satisfied, it would be nice if they considered something nonlinear as well.  This is especially true because they remove the regularization in their experiments anyhow and thus are at least somewhat happy to depart from their theory.

---

confidence: 5
contribution: 2
initial_rating: 1
presentation: 4
questions: n/a
soundness: 1
strengths: - The paper points out an important fact about the analysis of algorithms for "smoothed" spectral risk measures, which is that when the smoothing parameter $\nu = O(\epsilon)$, it should be included in the complexity guarantee, making "linearly convergent" methods sublinear for the original non-smooth problem.
- By using the experimental benchmark of [Mehta et al. (2024)](https://arxiv.org/pdf/2310.13863), the authors are able to perform head-to-head comparisons to  recent algorithms designed for spectral risk measures (Figure 3).  The empirical performance is competitive with, and often outperforms, other methods.
- The exposition is clear and well-written, although some sections are nearly identical to sections in the referenced work on spectral risk measure optimization with stochastic algorithms (Section 1, Section 2, Appx B).
summary: The paper presents SOREL, an optimization algorithm designed for spectral risk measure objectives, which can also be understood as  primal-dual min-max objectives for learning problems. In comparison to previous stochastic methods, the algorithm adds a proximal term to the dual update and claims a near-optimal convergence guarantee for the strongly convex-(nonstrongly) concave variant of this objective. The method employs an SVRG-style variance reduction similar to [Mehta et al. (2023)](https://proceedings.mlr.press/v206/mehta23b.html). While SOREL generally performs favorably in experiments, the theoretical convergence analysis is incorrect, making this work unpublishable in its current state.
weaknesses: Upon viewing Corollary 1, we see that the complexity contains the term $\log(\frac{\sqrt{n}G}{\mu^2\sqrt{\epsilon}})$, i.e. a logarithm is taken for a quantity that is *not* unitless. This does not pass a basic sanity check for convergence analyses in optimization theory: that the complexity does not depend on the units of the input. The error comes from the fact that the precision parameter $\delta_k$ in the inner loop is dependent on the strong convexity parameter $\mu$ (units of loss/inputs$^2$) when it should be in units of loss (based on Lemma 3). To see why this is invalid, recall that an optimization algorithm should not change its exact complexity when the inputs are reparametrized in different units. Consider the objective as a function of $\boldsymbol{w}/10$ instead of $\boldsymbol{w}$. Then, we have the Lipschitz and strong convexity constants change as $G \mapsto 10G$ and $\mu \mapsto 100\mu$. Then,
$$
\frac{\sqrt{n}(10G)}{(100\mu)^2\sqrt{\epsilon}} = \frac{1}{1000}\frac{\sqrt{n}G}{\mu^2\sqrt{\epsilon}}.
$$
This indicates that the number of iterations in the inner loop would change based on the scale, which is not a valid analysis. Please check the units before presenting the result.

Matters such as significance, novelty, and impact come second to the analysis, especially as the experimental results do not show a large practical improvement over the methods of [Mehta et al. (2023)](https://proceedings.mlr.press/v206/mehta23b.html) and [Mehta et al. (2024)](https://arxiv.org/pdf/2310.13863). Please correct the analysis.

---

confidence: 2
contribution: 3
initial_rating: 6
presentation: 3
questions: - You wrote in your abstract that there is a lack of guarantees for the original spectral risk. You said later in l.78 that SOREL is the first stochastic algorithm with convergence guarantees for the spectral risk minimisation problem. However, in Section 3.1 you directly consider a regularised version of the spectral risk with a strongly convex function (Eq. 1). Later, in Eq. (4) you added another regulariser in $\mathbf{\lambda}$ to avoid huge variations of $\lambda$ during the learning process. All of these reason suggests that you actually do NOT deal with the original spectral risk. Do I miss something?
- l. 268-269 : Those lines are unclear, I do not under whether there is a need to compute $\mathcal{O}(n\log(n))$ operations at line 5, within the loop of line 3 in Algorithm 1 or if this is the total number of computations. Can you explicitly write the time complexity of your algorithm and compare it with those proposed in literature? This would draw a useful theoretical link to understand why your algorithm is faster than others in Section 5.
- I am not sure to understand why Lemma 1 is stated here, as there is no proof sketch in the main document. What do I miss?
- In Theorem 1, can you please define $w^*$? I do not understand what it is.
- In Section 5, it seems that you focused exclusively on linear regression, would SOREL work with small neural nets with 1-2 hidden layers? Even if theoretical guarantees may not hold anymore, I am curious to know whether the methods remain tractable and efficient.
soundness: 3
strengths: - Sections 1, 2 and 3.1 are truly well-written, it has been a pleasure to read it
- The idea of using a primal dual method based on Lemma 1 of (Blondel et al. 2020) is clever and elegant
- Algorithmic performances are convincing and are accompanied by sound theory
- It is great to propose a mini batching version of SOREL
summary: This paper introduces SOREL, a stochastic gradient-based algorithm for spectral risk minimisation. After an introduction and an overview of this work, the authors carefully present the current drawbacks of current approaches of spectral risk minimisation alongside their new algorithm in Section 3. Section 4 provide convergence guarantees, and Section 5 gathers various experiments showing competitive performance and computational time with respect to existing methods.
weaknesses: See Questions part

Meta-review: An algorithm called SOREL is proposed for spectral risk minimization, with convergence analysis for convex losses together with a strongly convex regularization. Experiments was limited to linear models in the original submission, but additional ones with nonlinear models are conducted during the discussion. Revision should include the corrected proof and additional experiments.

Decision: Accept (Poster)