# Incremental Quasi-Newton Algorithms for Solving Nonconvex, Nonsmooth,
Finite-Sum Optimization Problems

Gulcin Dinc Yalcina and Frank E. Curtisb CONTACT G. Dinc Yalcin. Email:
<EMAIL_ADDRESS>ORCID ID: 0000-0001-7696-7507Frank E. Curtis ORCID ID:
0000-0001-7214-9187 aDepartment of Industrial Engineering, Eskisehir Technical
University, Eskisehir, 26555, Türkiye; a,bDepartment of Industrial and Systems
Engineering, Lehigh University, Bethlehem, PA 18015, USA

###### Abstract

Algorithms for solving nonconvex, nonsmooth, finite-sum optimization problems
are proposed and tested. In particular, the algorithms are proposed and tested
in the context of an optimization problem formulation arising in semi-
supervised machine learning. The common feature of all algorithms is that they
employ an incremental quasi-Newton (IQN) strategy, specifically an incremental
BFGS (IBFGS) strategy. One applies an IBFGS strategy to the problem directly,
whereas the others apply an IBFGS strategy to a difference-of-convex
reformulation, smoothed approximation, or (strongly) convex local
approximation. Experiments show that all IBFGS approaches fare well in
practice, and all outperform a state-of-the-art bundle method.

###### keywords:

nonconvex optimization, nonsmooth optimization, semi-supervised machine
learning, incremental quasi-Newton methods, smoothing

## 2010 MATHEMATICS SUBJECT CLASSIFICATIONS

49M37, 65K05, 90C26, 90C30, 90C53

##  1 Introduction

In this paper, we investigate the performance of algorithms for solving
nonconvex, nonsmooth, finite-sum minimization problems that arise in areas
such as statistical/machine learning. In particular, given a set of
$m\in\mathbb{N}$ functions that might be nonconvex and/or nonsmooth, say,
$f_{i}:{\cal R}\to\mathbb{R}$ for all $i\in\\{1,\dots,m\\}$ (where $\cal R$ is
some finite-dimensional real vector space), we are interested in algorithms
designed to minimize the potentially nonconvex and/or nonsmooth objective
$f:{\cal R}\to\mathbb{R}$ as in

$\min_{\omega\in\cal R}\ f(\omega),\ \ \text{where}\ \
f(\omega):=\tfrac{1}{m}\sum_{i=1}^{m}f_{i}(\omega).$ (1)

Specifically, to demonstrate the design and use of a suite of related
algorithms, we focus on a particular binary classification problem arising in
semi-supervised machine learning, which can be formulated as follows. Suppose
one is given a set of $p\in\mathbb{N}$ _labeled_ training samples
$\\{(x_{i},y_{i})\\}_{i=1}^{p}$ with $x_{i}\in\mathbb{R}^{n}$ and
$y_{i}\in\\{-1,1\\}$ for all $i\in\\{1,\dots,p\\}$ as well as a set of
$q\in\mathbb{N}$ _unlabeled_ training samples $\\{x_{i}\\}_{p+1}^{p+q}$ with
$x_{i}\in\mathbb{R}^{n}$ for all $i\in\\{p+1,\dots,p+q\\}$. The _transductive
support vector machine_ (TSVM) problem [6], similar to the classical support
vector machine (SVM) problem, aims to find a hyperplane that separates points
in the positively labeled class from those in the negatively labeled class
(labeled by $1$ and $-1$, respectively). An unconstrained optimization
formulation of the TSVM problem (see, e.g., [2, 37]), can be written in terms
of optimization variables $w\in\mathbb{R}^{n}$ and $b\in\mathbb{R}$ as the
problem to minimize $F:\mathbb{R}^{n}\times\mathbb{R}\to\mathbb{R}$, where

$F(w,b)=\tfrac{1}{2}\|w\|_{2}^{2}+C_{1}\sum_{i=1}^{p}\max\\{0,1-y_{i}(w^{T}x_{i}+b)\\}+C_{2}\sum_{i=p+1}^{p+q}\max\\{0,1-|w^{T}x_{i}+b|\\},$
(2)

with $C_{1}\in(0,\infty)$ and $C_{2}\in(0,\infty)$ being user-defined
parameters. In this definition of $F$, the first term is strongly convex, the
second term is convex and nonsmooth, and the third term is nonconvex and
nonsmooth; hence, $F$ is nonconvex and nonsmooth.

The algorithms that we propose and investigate are all built upon an
incremental quasi-Newton (IQN) framework, specifically the incremental version
of the Broyden-Fletcher-Goldfarb-Shanno (BFGS) [10, 27, 29, 48] strategy known
as IBFGS [42]. As discussed in Section 1.2, quasi-Newton methods (such as
BFGS) are designed to mimic the behavior of Newton’s method for optimization
through only the use of first-order derivative information, whereas
incremental quasi-Newton is designed to further exploit the finite-sum
structure of an objective function. Since our problem of interest, namely,
(2), is nonconvex and nonsmooth, the convergence guaranteees for the original
IBFGS algorithm [42] do not hold in our setting, although we still investigate
the performance of the approach applied directly to the problem to see if it
can still be efficient and effective in practice. We also propose and
investigate the performance of IBFGS strategies applied to solve a few
reformulations and approximations of problem (2), namely, one that uses a
difference-of-convex reformulation, one that uses a smoothed approximation,
and others that use (strongly) convex local approximations.

###  1.1 Contributions

Our main goal is to answer the question: Can IQN (such as IBFGS) strategies be
effective and efficient when solving the nonconvex and nonsmooth problem (2),
which would provide evidence that they should be considered viable options
more generally when solving problems of the form (1)? Our answer is in the
affirmative, as we show in numerical experiments that a straightforward
application of an IBFGS method, as well as application of all of our proposed
reformulations and approximations, are successful in our experiments, and all
outperform a state-of-the-art bundle method. Convergence guarantees for IBFGS
in our nonconvex and nonsmooth setting are elusive, as they are elusive for
numerous other algorithms for solving nonconvex and nonsmooth problems that
are applied regularly in practice to great effect, and as they are elusive
even for non-incremental quasi-Newton strategies in general in convex and
nonsmooth settings. Nevertheless, we show that our proposed strategies, which
involve either applying an IBFGS strategy directly or to a mollified
reformulation/approximation, are effective when solving challenging problems
motivated by an important real-world application.

###  1.2 Literature Review

Steepest descent methods (also known as gradient methods) for optimization
require only first-order derivative information and can achieve a local linear
rate of convergence in favorable settings. Newton methods, on the other hand,
use second-order derivative information to obtain a faster rate of local
convergence—namely, at least quadratic in favorable cases—but this comes at
the cost of additional computation to solve a large-scale linear system in
each iteration. Quasi-Newton methods offer a balance between per-iteration
computational cost and convergence rate by approximating second-order
derivatives using only first-order derivative information. The BFGS method is
one of the most popular quasi-Newton methods in practice [45]. Like other
quasi-Newton strategies, BFGS uses differences of gradients at consecutive
iterates to update iteratively an approximation of the (inverse) Hessian of
the objective. Convergence guarantees for the BFGS method have been proved for
convex smooth (see e.g. [13, 25]) and nonconvex smooth (see e.g. [40, 41])
objectives. BFGS methods have also been studied for the minimization of
nonsmooth functions (see, e.g., [22, 23]), although convergence guarantees for
nonsmooth minimization are very limited to certain convex problems (see [39]
for some example convergence guarantees for specific functions).

Steepest descent, Newton, and quasi-Newton methods have also been designed
specifically for the case of minimizing objectives defined by finite sums,
such as in (1). For example, incremental gradient methods use the gradient of
only a single component function in each iteration, cycling through the
component functions using a prescribed order or strategy. These methods can
significantly outperform non-incremental methods, although in theory they may
only achieve a sublinear rate of convergence [7]. The convergence properties
of such methods for convex minimization have been investigated in various
articles [8, 32, 50, 54]. For convex nonsmooth minimization, incremental
subgradient methods have been analyzed in, e.g., [36, 35, 38, 44, 47].

Beyond incremental (sub)gradient methods, [31] presents an incremental Newton
method for cases when the component functions are assumed to be strongly
convex. In [30], the convergence rate of incremental gradient and incremental
Newton methods are studied under constant and diminishing step size rules. For
the same reasons as in the non-incremental setting, IQN methods have also been
proposed for the incremental setting. After all, applying a method such as
BFGS when solving problem (1) requires ${\cal O}(mn+n^{2})$ operations per
iteration, which can be large when $m$ and/or $n$ is large. The IBFGS strategy
proposed in [42] has a cost of only ${\cal O}(n^{2})$ per iteration, and has
convergence guarantees when the objective is convex and smooth. It is this
algorithm that we use as a basis for the methods proposed in this paper.

Various other problems of the form (1), some related to (2), have been studied
in the literature, as well as algorithms for solving them; see, e.g., [37, 16,
5, 49, 20, 2]. We also refer the reader to review papers and book chapters,
such as [15, 51, 56]. This being said, to the best of our knowledge, there has
not been prior work on the variants of the IBFGS strategy that are proposed
and investigated in this paper. We also mention in passing that stochastic
(quasi-)Newton methods have been studied for solving machine learning problems
as well [11, 12, 46, 55].

###  1.3 Organization

In Section 2, we present reformulations and approximations of (2) to which our
proposed algorithms are to be applied. In Section 3, we provide background on
quasi-Newton methods, BFGS, and IBFGS, then propose variants of IBFGS for
solving the problems formulations that are presented in Section 2. Section 4
provides the results of numerical experiments using our proposed algorithms
and a state-of-the-art bundle method for the sake of comparison. Concluding
remarks are given in Section 5.

##  2 Reformulations and Approximations

As mentioned, the semi-supervised machine learning objective function (2)
involves three terms—one strongly convex, one convex and nonsmooth, and one
nonconvex and nonsmooth—leading overall to a function that is nonconvex and
nonsmooth. Like other such functions that appear in the literature, however,
this function has structure that can be exploited in the context of employing
optimization algorithms. In this section, we present a few potential
reformulations and approximations of this function that could be useful in
practice; indeed, our experiments confirm that they are useful.

###  2.1 Difference-of-Convex (DC) Function Reformulation

If one is willing to accept the nonsmoothness of (2), then the remaining
obstacle is nonconvexity. Fortunately, (2) has the property that it can be
expressed as a special type of nonconvex function, namely, a difference-of-
convex (DC) function, which is a type that is now well known to have
properties that can be exploited in the context of optimization [21].
Formally, let us introduce the following definition.

######  Definition 2.1.

(see, e.g., [34, Definition 2.1]) _A continuous function $f:\cal
R\to\mathbb{R}$ is a difference-of-convex (DC) function if and only if there
exist two convex functions, call them $g:\cal R\to\mathbb{R}$ and $h:\cal
R\to\mathbb{R}$, such that $f=g-h$._

Since the first two terms in (2) are convex, to show that $F$ is a DC
function, all that is needed is to express the third term as a DC function.
This can be done by recognizing that, for all $i\in\\{p+1,\dots,p+q\\}$, one
can write

$\displaystyle f_{i}(w,b):=$ $\displaystyle\ \max\\{0,1-|w^{T}x_{i}+b|\\}$
$\displaystyle=$ $\displaystyle\
\max\\{0,|w^{T}x_{i}+b|-1\\}-(|w^{T}x_{i}+b|-1)=:g_{i}(w,b)-h_{i}(w,b).$

This decomposition is illustrated in Figure 1. Using this observation, a DC
reformulation of (2) is given by
$F^{DC}:\mathbb{R}^{n}\times\mathbb{R}\to\mathbb{R}$ defined by

$F^{DC}(w,b)=\tfrac{1}{2}\|w\|_{2}^{2}+C_{1}\sum_{i=1}^{p}\max\\{0,1-y_{i}(w^{T}x_{i}+b)\\}\\\
+C_{2}\left(\sum_{i=p+1}^{p+q}\max\\{|w^{T}x_{i}+b|-1,0\\}-\sum_{i=p+1}^{p+q}(|w^{T}x_{i}+b|-1)\right).$
(3)

Figure 1: DC function formulation of $f_{i}(\cdot)=\max\\{0,1-|\cdot|\\}$ (on
the left) as the difference between $g_{i}(\cdot)=\max\\{0,|\cdot|-1\\}$ (in
the middle) and $h_{i}(\cdot)=|\cdot|-1$ (on the right).

###  2.2 Smooth Approximation

With the idea of opening the door for the use of the numerous algorithms that
have been designed for smooth optimization, one can consider avoiding the
nonsmoothness of (2) by considering a smooth approximation. In particular, one
can consider attempting to minimize (2) by solving (approximately) a sequence
of smooth approximations that tend to (2) in the limit. For this purpose, let
us introduce the following definition.

######  Definition 2.2.

(see, e.g., [17, Definition 1]) _A function $f^{S}:{\cal
R}\times(0,\infty)\to\mathbb{R}$ is a smoothing function of a continuous
function $f:\cal R\to\mathbb{R}$ if and only if $f^{S}(\cdot,\mu)$ is
continuously differentiable for all $\mu\in(0,\infty)$ and, for any
$\omega\in\cal R$, one has_

$\lim_{\bar{\omega}\to\omega,\mu\searrow 0}f^{S}(\bar{\omega},\mu)=f(\omega).$

Fortunately, the nonsmoothness of (2) can be smoothed in a straightforward
manner to obtain a smoothing function. Following [17, 18], for the univariate
function $(t)_{+}=\max\\{0,t\\}$, a smoothing function
$\phi:\mathbb{R}\times\mathbb{R}\to\mathbb{R}$ is given by either

$\phi(t,\mu)=\begin{cases}(t)_{+}&\text{if}\quad|t|\geq\tfrac{\mu}{2}\\\
\tfrac{t^{2}}{2\mu}+\tfrac{t}{2}+\tfrac{\mu}{8}&\text{if}\quad|t|<\tfrac{\mu}{2}\end{cases}\
\ \text{or}\ \ \phi(t,\mu)=\tfrac{1}{2}(t+\sqrt{t^{2}+4\mu^{2}}).$ (4)

Correspondingly, using the fact that $|t|=(t)_{+}+(-t)_{+}$, a smoothing
function $\psi:\mathbb{R}\times\mathbb{R}\to\mathbb{R}$ for the univariate
absolute value function is given by either

$\psi(t,\mu)=\begin{cases}|t|&\text{if}\quad|t|\geq\tfrac{\mu}{2}\\\
\tfrac{t^{2}}{\mu}+\tfrac{\mu}{4}&\text{if}\quad|t|<\tfrac{\mu}{2}\end{cases}\
\ \text{or}\ \ \psi(t,\mu)=\sqrt{t^{2}+4\mu^{2}}.$ (5)

Other choices are possible as well, each with different theoretical and/or
practical advantages and disadvantages. For example, the former options above
have the advantage that they match $(t)_{+}$ and $|t|$, respectively, when
$|t|\geq\tfrac{\mu}{2}$, which is not the case for the latter options.
However, the latter options have the advantage of being twice continuously
differentiable, which is not the case for the former options.

One can now obtain a smooth approximation of (2) using smoothing functions.
First, for all $i\in\\{1,\dots,p\\}$, suppose a smoothing function $f_{i}^{S}$
of $f_{i}$ is given by $f_{i}^{S}(w,b,\mu)=\phi(t_{i},\mu)$, where
$t_{i}\equiv 1-y_{i}(w^{T}x_{i}+b)$. Second, for all
$i\in\\{p+1,\dots,p+q\\}$, suppose a smoothing function $f_{i}^{S}$ of $f_{i}$
is given by $f_{i}^{S}(w,b,\mu)=\phi(t_{i},\mu)$, where $t_{i}\equiv
1-\psi(w^{T}x_{i}+b,\mu)$. Then, a smoothing function $F^{S}$ for $F$ in (2)
is given by

$\displaystyle F^{S}(w,b,\mu)=\tfrac{1}{2}\|w\|_{2}^{2}$
$\displaystyle+C_{1}\sum_{i=1}^{p}\phi(1-y_{i}(w^{T}x_{i}+b),\mu)$ (6)
$\displaystyle+C_{2}\sum_{i=p+1}^{p+q}\phi(1-\psi(w^{T}x_{i}+b,\mu),\mu);$

see Figure 2 for illustration. Here, we are making use of the fact that a
smooth composite of smoothing functions is a smoothing function; see, e.g.,
[17, Theorem 1].

Figure 2: Illustration of $\phi(1-(\cdot),\mu)$ (left) and
$\phi(1-\psi(\cdot,\mu),\mu)$ (right) with $\mu=0.1$ using the former choices
for the smoothing functions in (4) and (5).

###  2.3 Smooth Convex Local Approximation

The smoothing function (6) has the benefit of being continuously
differentiable. However, it is nonconvex, which presents certain challenges
when one aims to employ optimization techniques such as a quasi-Newton
strategy. To address these challenges, one can consider—about some base point,
call it $(\bar{w},\bar{b})\in\mathbb{R}^{n}\times\mathbb{R}$—a convex “local”
approximation of $F^{S}$. Thus, the question becomes: Can $F^{S}$ be
“convexified” locally in a straightforward manner? The answer is yes, since
$F^{S}$ happens to be weakly convex.

######  Definition 2.3.

(see, e.g., [53, Definition 4.1 and Proposition 4.3]) _A function $f:{\cal
R}\to\mathbb{R}$ is $\rho$-convex if and only if there exists a convex
function $g:{\cal R}\to\mathbb{R}$ such that
$f(\omega)=g(\omega)+\rho\|x\|_{2}^{2}$. A function that is $\rho$-convex for
$\rho\in(0,\infty)$ is said to be strongly convex. A function that is
$\rho$-convex for $\rho\in(-\infty,0)$ is said to be weakly convex._

To see that $F^{S}$ is weakly convex, the following proposition is useful.

######  Proposition 2.4.

(see, e.g., [24, Section 3.1]) If $c:{\cal R}\to\mathbb{R}$ is continuously
differentiable with its gradient function being Lipschitz with constant
$L_{\nabla c}\in(0,\infty)$, and if $h:\mathbb{R}\to\mathbb{R}$ is convex and
Lipschitz with constant $L_{h}\in(0,\infty)$, then the composite function
$h\circ c:{\cal R}\to\mathbb{R}$ is $\rho$-convex with $\rho=-L_{h}L_{\nabla
c}\in(-\infty,0)$.

The first two terms for $F^{S}$ in (6) are convex. Hence, to see that $F^{S}$
is weakly convex, it follows from Proposition 2.4 that it suffices to notice
that the final term is a composite function following the conditions of the
proposition. For example, let us consider the former options in (4) and (5).
One finds for any $\mu\in(0,\infty)$ that $\phi(\cdot,\mu)$ is convex and
Lipschitz with constant 1. In addition, one finds for any $\mu\in(0,\infty)$
and $i\in\\{p+1,\dots,p+q\\}$ that $1-\psi(w^{T}x_{i}+b,\mu)$ is continuously
differentiable with a gradient that is Lipschitz continuous with constant
$2(\|x_{i}\|_{2}^{2}+1)/\mu$. Hence, given any smoothing parameter
$\mu\in(0,\infty)$ and base point $(\bar{w},\bar{b})$, a smooth and convex
local approximation of $F^{S}$ is given by
$F^{C}(\cdot,\cdot,\mu,\bar{w},\bar{b}):\mathbb{R}^{n}\times\mathbb{R}\to\mathbb{R}$
defined by

$F^{C}(w,b,\mu,\bar{w},\bar{b})=\tfrac{1}{2}\|w\|_{2}^{2}+C_{1}\sum_{i=1}^{p}\phi_{i}(1-y_{i}(w^{T}x_{i}+b),\mu)\\\
+C_{2}\sum_{i=p+1}^{p+q}\left(\phi(1-\psi(w^{T}x_{i}+b,\mu),\mu)+\rho_{i}\left\|\begin{bmatrix}w\\\
b\end{bmatrix}-\begin{bmatrix}\bar{w}\\\
\bar{b}\end{bmatrix}\right\|_{2}^{2}\right),$ (7)

where for all $i\in\\{p+1,\dots,p+q\\}$ one has chosen $\rho_{i}\in(0,\infty)$
with $\rho_{i}\geq 2(\|x_{i}\|_{2}^{2}+1)/\mu$.

Notice that it has been critical in this discussion to consider the smooth
approximation $F^{S}(\cdot,\mu)$ for some $\mu\in(0,\infty)$, rather than $F$
in (2) directly. After all, the function $F$ in (2) is not weakly convex.
Therefore, it should not be surprising that the lower bounds for
$\\{\rho_{i}\\}_{i=p+1}^{p+q}$ tend to infinity as $\mu\searrow 0$, i.e., as
$F^{S}$ tends to $F$.

###  2.4 Smooth Strongly Convex Local Approximation

The algorithms that we propose and test in the following sections are built on
the fact that the objective functions under consideration consist of a finite
sum of functions. For example, at face value, such an algorithm could consider
$\tfrac{1}{2}\|w\|_{2}^{2}$ as one function and the remaining $p+q$ terms as
separate functions, meaning that one can consider there being $p+q+1$ separate
terms overall. However, one disadvantage of this breakdown is that not all of
these terms would be strongly convex. This is a disadvantage since, as we
shall mention in the next section, quasi-Newton methods can avoid the need to
_damp_ or _skip_ updates when a function is strongly convex. Hence, it is
worthwhile also to explore an alternative to (7) in which all terms are
strongly convex.

One way to derive such an alternative is to simply add a strongly convex
quadratic function to all terms. However, this turns out not to be necessary.
Instead, one can consider _distributing_ the $\tfrac{1}{2}\|w\|_{2}^{2}$ term
into the remaining $p+q$ terms, which has a similar effect, at least in terms
of the $w$ variables. Merely adding and also distributing a term of the form
$\tfrac{1}{2}\beta\|b-\bar{b}\|_{2}^{2}$, one finds that the following only
differs from (7) by $\tfrac{1}{2}\beta\|b-\bar{b}\|_{2}^{2}$ and has the
advantage that each of the $p+q$ terms is strongly convex:

$F^{SC}(w,b,\mu,\bar{w},\bar{b})=C_{1}\sum_{i=1}^{p}(\phi(1-y_{i}(w^{T}x_{i}+b),\mu)+\tfrac{1}{4pC_{1}}\|w\|_{2}^{2}+\tfrac{\beta}{2pC_{1}}\|b-\bar{b}\|_{2}^{2})\\\
+C_{2}\sum_{i=p+1}^{p+q}(\phi(1-\psi(w^{T}x_{i}+b,\mu),\mu)+\tfrac{\rho_{i}}{2}(\|w-\bar{w}\|_{2}^{2}+\|b-\bar{b}\|_{2}^{2})+\tfrac{1}{4qC_{2}}\|w\|_{2}^{2}),$
(8)

where $\beta\in(0,\infty)$ and
$\rho_{i}\in(2(\|x_{i}\|_{2}^{2}+1)/\mu,\infty)$ for all
$i\in\\{p+1,\dots,p+q\\}$. (Again, we are assuming that the former options in
(4) and (5) are being used. The use of alternatives would only potentially
affect the allowable ranges for $\\{\rho_{i}\\}_{i=p+1}^{p+q}$.)

Formally, the theoretical benefits of having terms that are strongly convex
are that they satisfy certain inequalities that, as we shall see in the next
section, are required by the known theoretical convergence guarantees of
IBFGS.

######  Lemma 2.5.

If a function $f:{\cal R}\to\mathbb{R}$ is continuously differentiable over a
sublevel set ${\cal L}_{\gamma}(f):=\\{x\in{\cal R}:f(x)\leq\gamma\\}$ for
some $\gamma\in\mathbb{R}$, its gradient function is Lipschitz continuous with
constant $L_{\nabla f}$ over ${\cal L}_{\gamma}(f)$, and $f$ is strongly
convex over ${\cal L}_{\gamma}(f)$, then there exists $\mu_{\nabla
f}\in(0,\infty)$ such that, for all $(\omega,\bar{\omega})\in{\cal
L}_{\gamma}(f)\times{\cal L}_{\gamma}(f)$, one has

$\mu_{\nabla f}\|\omega-\bar{\omega}\|_{2}^{2}\leq(\nabla f(\omega)-\nabla
f(\bar{\omega}))^{T}(\omega-\bar{\omega})\leq L_{\nabla
f}\|\omega-\bar{\omega}\|_{2}^{2}.$

###### Proof.

The former inequality follows from [53, Proposition 4.10]. The latter follows
from the Lipschitz continuity of the gradient of $f$ over ${\cal
L}_{\gamma}(f)$. ∎

##  3 (Incremental) Quasi-Newton Methods

Suppose that one aims to minimize a twice-continuously differentiable
objective function. Quasi-Newton methods attempt to mimic the behavior of
Newton’s method for the minimization of such a function using only first-order
derivative information. This is accomplished by maintaining a matrix
representing an approximation of the Hessian of the objective. This matrix is
updated iteratively using iterate and gradient displacements during the
minimization algorithm, and in favorable situations the minimization algorithm
can attain a superlinear rate of convergence [9].

Mathematically, a quasi-Newton strategy operates as follows. For simplicity of
notation, let us denote the objective function as
$f:\mathbb{R}^{n+1}\to\mathbb{R}$ (i.e., let us update our notation so that
${\cal R}\equiv\mathbb{R}^{n+1}$ where, for our purposes later on, we replace
$(w,b)\in\mathbb{R}^{n}\times\mathbb{R}$ in a natural manner with
$\omega\equiv(w^{T},b^{T})^{T}\in\mathbb{R}^{n+1}$) and suppose that the
minimization algorithm has reached iteration $k\in\mathbb{N}$ and that the
current approximation is the symmetric positive definite matrix
$B_{k}\in\mathbb{R}^{n+1}\times\mathbb{R}^{n+1}$. It is assumed that the
minimization algorithm computes a search direction $p_{k}$ as the minimizer of
the local quadratic approximation of $f$ at $\omega_{k}\in\mathbb{R}^{n+1}$
given by $\nabla f(\omega_{k})^{T}p+\tfrac{1}{2}p^{T}B_{k}p$, computes a step
size $\alpha_{k}\in(0,\infty)$, then sets the next iterate as
$\omega_{k+1}\leftarrow\omega_{k}-\alpha_{k}B_{k}^{-1}\nabla f(\omega_{k})$,
as in any Newton-like scheme. The signifying feature of a quasi-Newton
strategy is the manner in which the subsequent approximation matrix
$B_{k+1}\in\mathbb{R}^{n+1}\times\mathbb{R}^{n+1}$ is chosen. In particular,
in a quasi-Newton strategy, one defines the iterate and gradient displacements

$s_{k}=\omega_{k+1}-\omega_{k}\ \ \text{and}\ \ y_{k}=\nabla
f(\omega_{k+1})-\nabla f(\omega_{k}),$ (9)

respectively, then chooses $B_{k+1}$ satisfying the so-called secant equation
$B_{k+1}s_{k}=y_{k}$.

The secant equation does not determine $B_{k+1}$ uniquely; there are various
choices that satisfy it. The popular and effective BFGS method makes the
selection of the matrix satisfying the secant equation whose inverse minimizes
the difference from the inverse of the prior approximation. In particular, the
BFGS method chooses $B_{k+1}$ as the unique solution to the optimization
problem

$\min_{B^{-1}}\|B^{-1}-B_{k}^{-1}\|\ \ \text{s.t.}\ \ B^{-1}=B^{-T}\ \
\text{and}\ \ B^{-1}y_{k}=s_{k},$

where the norm in the objective of this problem is a particular weighted
Frobenius norm chosen such that a rank-two change in $B_{k}^{-1}$ yields
$B_{k+1}^{-1}$ [45]. Applying the Sherman-Morrison-Woodbury formula to the
formula for $B_{k+1}^{-1}$ in terms of $B_{k}^{-1}$, one obtains the following
formula for $B_{k+1}$ in terms of $B_{k}$:

$B_{k+1}=B_{k}-\frac{B_{k}s_{k}s_{k}^{T}B_{k}}{s_{k}^{T}B_{k}s_{k}}+\frac{y_{k}y_{k}^{T}}{y_{k}^{T}s_{k}}.$

For future reference, we remark that if the pair $(s_{k},y_{k})$ satisfies
$s_{k}^{T}y_{k}>0$ for all $k\in\mathbb{N}$, then as long as the initial
approximation $B_{0}$ is symmetric and positive definite (e.g., $B_{0}=I$),
then all elements of the sequence $\\{B_{k}\\}$ will be symmetric and positive
definite as well. Moreover, $s_{k}^{T}y_{k}>0$ is guaranteed to hold for all
$k\in\mathbb{N}$ if the gradients employed in the definition of $y_{k}$ are
computed from a strongly convex function [45].

If the objective function $f$ is not smooth, but still locally Lipschitz, then
while standard convergence guarantees for quasi-Newton methods no longer
apply, it has been shown to be effective in practice to employ quasi-Newton
methods with generalized gradients in place of gradients. Following the
developments of Clarke [19], the set of generalized gradients for a locally
Lipschitz function $f$ at $\omega$ is given by

$\eth f(\omega)=\text{conv}\left(\left\\{\lim_{\bar{\omega}\to\omega,\
\bar{\omega}\in{\cal D}_{f}}\nabla f(\omega_{k})\right\\}\right),$

where ${\cal D}_{f}\subseteq\mathbb{R}^{n+1}$ is the set of points at which
$f$ is differentiable. For many problems, at any $\omega$, one is able to
compute an (arbitrary) element from $\eth f(\omega)$, meaning that in a
practical application of a quasi-Newton method one can replace $y_{k}$ from
(9) with

$y_{k}=g_{k+1}-g_{k},\ \ \text{where}\ \ (g_{k+1},g_{k})\in\eth
f(\omega_{k+1})\times\eth f(\omega_{k}).$

###  3.1 I-BFGS

Incremental algorithms, as described in Section 1.2, attempt to exploit the
finite-sum structure of an objective function by using information from only a
single component function to compute a step during each iteration, cycling
through the component functions using a prescribed order or strategy.
Following this idea, for an incremental quasi-Newton method for the
minimization of a finite-sum objective function of the form (1), an idea is to
store a Hessian approximation matrix for each component function individually.
Then, since new information is obtained for only a single component function
in each iteration, the Hessian approximation for that component function can
be updated, e.g., using BFGS, while all others are left the same [42].

An I-BFGS strategy of this type can be described as follows. For all
$k\in\mathbb{N}$ and $i\in\\{1,\dots,m\\}$, let $z_{k,i}\in\mathbb{R}^{n+1}$
denote the latest point at which a generalized gradient has been computed for
component function $i$ up to iteration $k$, let $v_{k,i}\in\eth
f_{i}(z_{k,i})$ denote this latest generalized gradient for component function
$i$ up to iteration $k$, and let $B_{k,i}\in\mathbb{R}^{n+1\times n+1}$ denote
the Hessian approximation for component function $i$ in iteration $k$. With
respect to these quantities, the algorithm initializes
$z_{0,i}\leftarrow\omega_{0}$, $v_{0,i}\in\eth f_{i}(z_{0,i})$, and
$B_{0,i}\succ 0$ and for all $i\in\\{1,\dots,m\\}$. During iteration $k$, a
new iterate is computed using information available at the start of iteration
$k$, then the algorithm chooses an index $i_{k}\in\\{1,\dots,m\\}$
corresponding to which the Hessian approximation will be updated. In
particular, the algorithm sets $z_{k+1,i_{k}}\leftarrow\omega_{k+1}$,
$v_{k+1,i_{k}}\in\eth f_{i_{k}}(z_{k+1,i_{k}})$, $s_{k,i_{k}}\leftarrow
z_{k+1,i_{k}}-z_{k,i_{k}}$, and $y_{k,i_{k}}\leftarrow
v_{k+1,i_{k}}-v_{k,i_{k}}$, then sets

$B_{k+1,i_{k}}\leftarrow
B_{k,i_{k}}-\frac{B_{k,i_{k}}s_{k,i_{k}}s_{k,i_{k}}^{T}B_{k,i_{k}}}{s_{k,i_{k}}^{T}B_{k,i_{k}}s_{k,i_{k}}}+\frac{y_{k,i_{k}}y_{k,i_{k}}^{T}}{y_{k,i_{k}}^{T}s_{k,i_{k}}},$
(10)

whereas, for all other indices, the points, generalized gradients, and Hessian
approximations are unchanged, and the iterate and gradient displacements are
set to zero.

In a straightforward application of this idea, in iteration $k\in\mathbb{N}$,
one is left with the desire to compute the new iterate $\omega_{k+1}$ as the
minimizer of $q_{k}:\mathbb{R}^{n+1}\to\mathbb{R}$ defined as

$q_{k}(\omega)=\frac{1}{m}\sum_{i=1}^{m}\left(f_{i}(z_{k,i})+v_{k,i}^{T}(w-z_{k,i})+\frac{1}{2}(w-z_{k,i})^{T}B_{k,i}(w-z_{k,i})\right),$

the solution of which can be expressed as

$\omega_{k+1}=\left(\sum_{i=1}^{m}B_{k,i}\right)^{-1}\left(\sum_{i=1}^{m}B_{k,i}z_{k,i}-\sum_{i=1}^{m}v_{k,i}\right).$

However, forming $\sum_{i=1}^{m}B_{k,i}$ by summing the elements in
$\\{B_{k,i}\\}_{i=1}^{m}$, then using the resulting sum to compute
$\omega_{k+1}$ could be prohibitively expensive. Instead, in a more
computationally efficient approach, the matrix
$(\tilde{B}_{k})^{-1}:=(\sum_{i=1}^{m}B_{k,i})^{-1}$ can be updated directly.
In particular, following [42], one can set

$(\tilde{B}_{k+1})^{-1}\leftarrow
U_{k}+\frac{U_{k}(B_{k,i_{k}}s_{k,i_{k}})(B_{k,i_{k}}s_{k,i_{k}})^{T}U_{k}}{s_{k,i_{k}}^{T}B_{k,i_{k}}s_{k,i_{k}}-(B_{k,i_{k}}s_{k,i_{k}})^{T}U_{k}(B_{k,i_{k}}s_{k,i_{k}})},$
(11)

where

$U_{k}\leftarrow(\tilde{B}_{k})^{-1}+\frac{(\tilde{B}_{k})^{-1}y_{k,i_{k}}y_{k,i_{k}}^{T}(\tilde{B}_{k})^{-1}}{y_{k,i_{k}}^{T}s_{k,i_{k}}+y_{k,i_{k}}^{T}(\tilde{B}_{k})^{-1}y_{k,i_{k}}}.$
(12)

We state a complete algorithm, which follows that in [42], as Algorithm 1. We
remark that this algorithm employs unit step sizes, although in an
implementation one should consider potentially taking shorter step sizes to
ensure convergence or at least good practical behavior, as we do in our
numerical experiments. We also note that Algorithm 1 skips the BFGS update in
iteration $k\in\mathbb{N}$ unless the iterate and generalized gradient
displacements are sufficiently large and a curvature condition holds. This is
a standard safeguard that is needed in nonconvex settings.

Algorithm 1 I-BFGS Method

1:$c\in(0,\infty)$

2:choose $\omega_{0}\in\mathbb{R}^{n+1}$.

3:set $z_{0,i}\leftarrow\omega_{0}$, $v_{0,i}\in\eth f_{i}(z_{0,i})$, and
$B_{0,i}\leftarrow I$ for all $i\in\\{1,\dots,m\\}$

4:set $(\tilde{B}_{0})^{-1}\leftarrow(\sum_{i=1}^{m}B_{0,i})^{-1}$,
$u_{0}\leftarrow\sum_{i=1}^{m}B_{0,i}\omega_{0}$, and
$g_{0}=\sum_{i=1}^{m}v_{0,i}$

5:for $k\in\mathbb{N}$ do

6: set $\omega_{k+1}\leftarrow(\tilde{B}_{k})^{-1}(u_{k}-g_{k})$

7: choose $i_{k}\in\\{1,\dots,m\\}$

8: set $z_{k+1,i_{k}}\leftarrow\omega_{k+1}$ and $z_{k+1,i}\leftarrow z_{k,i}$
for all $i\in\\{1,\dots,m\\}\setminus\\{i_{k}\\}$

9: set $v_{k+1,i_{k}}\in\eth f_{i_{k}}(z_{k+1,i_{k}})$ and
$v_{k+1,i}\leftarrow v_{k,i}$ for all
$i\in\\{1,\dots,m\\}\setminus\\{i_{k}\\}$

10: set $s_{k,i}\leftarrow z_{k+1,i}-z_{k,i}$ for all $i\in\\{1,\dots,m\\}$

11: set $y_{k,i}\leftarrow v_{k+1,i}-v_{k,i}$ for all $i\in\\{1,\dots,m\\}$

12: if $s_{k,i_{k}}^{T}y_{k,i_{k}}>c\|s_{k,i_{k}}\|_{2}\|y_{k,i_{k}}\|_{2}$,
$\|s_{k,i_{k}}\|_{2}>c$, and $\|y_{k,i_{k}}\|_{2}>c$ then

13: set $B_{k+1,i_{k}}$ by (10) and $B_{k+1,i}\leftarrow B_{k,i}$ for all
$i\in\\{1,\dots,m\\}\setminus\\{i_{k}\\}$

14: set $(\tilde{B}_{k+1})^{-1}$ by (11)–(12).

15: else

16: set $B_{k+1,i}\leftarrow B_{k,i}$ for all $i\in\\{1,\dots,m\\}$

17: set $(\tilde{B}_{k+1})^{-1}\leftarrow(\tilde{B}_{k})^{-1}$

18: end if

19: set $u_{k+1}\leftarrow
u_{k}+(B_{k+1,i_{k}}\omega_{k+1}-B_{k,i_{k}}z_{k,i_{k}})$

20: set $g_{k+1}\leftarrow g_{k}+v_{k+1,i_{k}}-v_{k,i_{k}}$

21:end for

In the following subsections, we state our proposed modifications of I-BFGS
that exploit the problem reformulations and approximations described in
Section 2.

###  3.2 I-BFGS-DC for DC Function Reformulation

As mentioned in Section 2.1, the structure of a DC function can be exploited
in the context of an optimization algorithm. For example, in the context of an
I-BFGS strategy, one can exploit the properties of a DC function in an attempt
to avoid the skipping of updates. Specifically, suppose that $f_{i}$ is a DC
function that can be expressed as $f_{i}=g_{i}-h_{i}$, where $g_{i}$ and
$h_{i}$ are convex functions. If, corresponding to a pair of points
$(\omega,\bar{\omega})$, one finds that the displacements
$s=\bar{\omega}-\omega$ and $y=v_{\bar{\omega}}-v_{\omega}$, where
$(v_{\omega},v_{\bar{\omega}})\in\eth f_{i}(\omega)\times\eth
f_{i}(\bar{\omega})$, yield $s^{T}y<0$, then a rule such as that employed in
Algorithm 1 would skip the BFGS update. This might be avoided by replacing $y$
by an approximate generalized gradient displacement, where the approximation
exploits knowledge of the DC function structure of $f_{i}$.

The strategy that we propose is the following. Suppose with respect to
$\omega$ one has

$v_{g_{i},\omega}\in\eth g_{i}(\omega),\ \ v_{h_{i},\omega}\in\eth
h_{i}(\omega),\ \ \text{and}\ \ v_{\omega}=v_{g_{i},\omega}-v_{h_{i},\omega}.$

Then, instead of evaluating $v_{\bar{\omega}}$ in a similar manner using a
generalized gradient of $h_{i}$ at $\bar{\omega}$ (i.e., a subgradient, since
$h_{i}$ is convex), one can instead employ a generalized gradient of the
affine underestimator of $h_{i}$ defined by
$h_{i}(\omega)+v_{h_{i},\omega}^{T}(\cdot)$, leading to

$v_{g_{i},\bar{\omega}}\in\eth g_{i}(\bar{\omega}),\ \
v_{h_{i},\bar{\omega}}\leftarrow v_{h_{i},\omega},\ \ \text{and}\ \
v_{\bar{\omega}}=v_{g_{i},\bar{\omega}}-v_{h_{i},\bar{\omega}},$

which in turn leads to the approximate displacement

$y=v_{\bar{\omega}}-v_{\omega}=v_{g_{i},\bar{\omega}}-v_{g_{i},\omega}.$

Since $g_{i}$ is convex, this means that $s^{T}y\geq 0$, and it is more likely
that $s^{T}y>0$.

The algorithm that we refer to as I-BFGS-DC is identical to I-BFGS (Algorithm
1), except that each instance of $\eth f_{i}(\cdot)$ is replaced by $\eth
g_{i}(\cdot)-\eth h_{i}(\cdot)$ (where “$-$” denotes a Minkowski difference)
in all places besides the computation of $v_{k+1,i_{k}}$ in Step 9, where it
is replaced by the subroutine stated in Algorithm 2. We include the if
condition in the subroutine since this alternative computation of
$v_{k+1,i_{k}}$ is not needed when $f_{i}$ is known to be convex. For example,
in I-BFGS-DC employed to minimize (2), this modification would only be
employed if $i_{k}$ corresponds to one of the final $q$ terms; it is not
needed for any of the other terms since they are convex. In any case, after
the new BFGS approximation is set, the subroutine states that the value of
$v_{k+1,i_{k}}$ is “reset” (if necessary) to an “unmodified” value for use in
later iterations.

Algorithm 2 I-BFGS-DC modified of computation of $v_{k+1,i_{k}}$

if $f_{i_{k}}$ is nonconvex then

set $v_{k+1,i_{k}}\in\eth g_{i_{k}}(z_{k+1,i_{k}})-\eth
h_{i_{k}}(z_{k,i_{k}})$

else

set $v_{k+1,i_{k}}\in\eth g_{i_{k}}(z_{k+1,i_{k}})-\eth
h_{i_{k}}(z_{k+1,i_{k}})$

end if

[… all proceeds as in I-BFGS until after line 17, then…]

if $f_{i_{k}}$ is nonconvex then

set $v_{k+1,i_{k}}\in\eth g_{i_{k}}(z_{k+1,i_{k}})-\eth
h_{i_{k}}(z_{k+1,i_{k}})$

end if

###  3.3 I-BFGS-S for Smooth Approximation

I-BFGS can be applied to a smooth approximation of (1) in a straightforward
manner. In particular, for all $i\in\\{1,\dots,m\\}$, letting
$f_{i}^{S}(\cdot,\mu_{i}):\mathbb{R}^{n+1}\to\mathbb{R}$ denote a smoothing
function approximation of $f_{i}$ given a smoothing parameter
$\mu_{i}\in(0,\infty)$, problem (1) (with ${\cal R}=\mathbb{R}^{n+1}$, as has
been used in this section) is approximated by

$\min_{\omega\in\mathbb{R}^{n+1}}\tfrac{1}{m}\sum_{i=1}^{m}f_{i}^{S}(\omega,\mu_{i}).$
(13)

The algorithm that we refer to as I-BFGS-S is identical to I-BFGS (Algorithm
1) except that the algorithm initializes $\mu_{0,i}\leftarrow\mu$ for all
$i\in\\{1,\dots,m\\}$ for some $\mu\in(0,\infty)$, each instance of an element
of $\eth f_{i}(\cdot)$ is replaced by $\nabla f_{i}(\cdot,\mu_{k,i})$, and the
smoothing parameter update in Algorithm 3 is added immediately after Step 10;
see [17].

Algorithm 3 I-BFGS-S addition for updating smoothing parameters

$\kappa\in(0,\infty)$, $\sigma\in(0,1)$

set $\mu_{k+1,i}\leftarrow\mu_{k,i}$ for all
$i\in\\{1,\dots,m\\}\setminus\\{i_{k}\\}$

if $\|v_{k+1,i_{k}}\|_{2}\geq\kappa\mu_{k,i_{k}}$ then

set $\mu_{k+1,i_{k}}\leftarrow\mu_{k,i_{k}}$

else

set $\mu_{k+1,i_{k}}\leftarrow\sigma\mu_{k,i_{k}}$ and replace
$v_{k+1,i_{k}}\leftarrow\nabla f_{i_{k}}^{S}(z_{k+1,i_{k}},\mu_{k+1,i_{k}})$

end if

We remark that since the gradients computed in I-BFGS-S depend on the
smoothing parameters, it is possible that a gradient displacement—employed in
a BFGS update—involves a difference of gradients computed using two different
values of the smoothing parameter. We did not find this to be an issue in our
experiments, but it is something of which one should be aware in the use of
I-BFGS-S.

###  3.4 I-BFGS-C for Smooth Convex Local Approximation

Suppose that, for all $i\in\\{1,\dots,m\\}$, the function
$f_{i}^{S}(\cdot,\mu_{i})$ is a continuously differentiable and
$\rho_{i}$-convex approximation of $f_{i}$. For example, for the function
$F^{S}$ in (6), the first term is continuously differentiable and
$1/2$-convex, the next $p$ terms are continuously differentiable and
$0$-convex, and the remaining $q$ terms are continuously differentiable and
$\rho_{i}$-convex for sufficiently negative $\rho_{i}$, as shown in Section
2.3. Then, in the context of an I-BFGS-type algorithm, rather than compute the
gradient displacement for the function with index $i_{k}$ as in I-BFGS-S, one
can attempt to avoid the need to skip an update by computing the displacement
according to a smooth convex local approximation of $f_{i_{k}}$. This can be
done following the discussion in Section 2.3 with the “base point” always set
as the last iterate at which a gradient was computed.

The algorithm that we refer to as I-BFGS-C is identical to I-BFGS (Algorithm
1) except that the algorithm initializes $\mu_{0,i}\leftarrow\mu$ for all
$i\in\\{1,\dots,m\\}$ for some $\mu\in(0,\infty)$, each instance of an element
of $\eth f_{i}(\cdot)$ is replaced by $\nabla f_{i}(\cdot,\mu_{k,i})$, and
Step 11 is replaced by the subroutine stated as Algorithm 4. This subroutine
updates the smoothing parameters using the same strategy as in I-BFGS-S, and
modifies the update for the gradient displacement if the smooth and
$\rho_{i_{k}}$-convex approximation of $f_{i_{k}}$ is nonconvex.

Algorithm 4 I-BFGS-C additions

1:$\kappa\in(0,\infty)$, $\sigma\in(0,1)$

2:set $\mu_{k+1,i}\leftarrow\mu_{k,i}$ for all
$i\in\\{1,\dots,m\\}\setminus\\{i_{k}\\}$

3:set $y_{k,i}\leftarrow v_{k+1,i}-v_{k,i}$ for all
$i\in\\{1,\dots,m\\}\setminus\\{i_{k}\\}$

4:if $\|v_{k+1,i_{k}}\|_{2}\geq\kappa\mu_{k,i_{k}}$ then

5: set $\mu_{k+1,i_{k}}\leftarrow\mu_{k,i_{k}}$

6:else

7: set $\mu_{k+1,i_{k}}\leftarrow\sigma\mu_{k,i_{k}}$ and replace
$v_{k+1,i_{k}}\leftarrow\nabla f_{i_{k}}^{S}(z_{k+1,i_{k}},\mu_{k+1,i_{k}})$

8:end if

9:if $f_{i_{k}}^{S}(\cdot,\mu_{k+1,i_{k}})$ is convex then

10: set $y_{k,i_{k}}\leftarrow v_{k+1,i_{k}}-v_{k,i_{k}}$

11:else

12: set $\rho$ such that
$f_{i_{k}}^{S}(\cdot,\mu_{k+1,i_{k}})+\rho\|\cdot\|_{2}^{2}$ is convex

13: set $y_{k,i_{k}}\leftarrow
v_{k+1,i_{k}}-v_{k,i_{k}}+2\rho(z_{k+1,i_{k}}-z_{k,i_{k}})$

14:end if

###  3.5 I-BFGS-SC for Smooth Strongly Convex Local Approximation

The algorithm that we refer to as I-BFGS-SC is identical to I-BFGS-C except
that, in line 12 of Algorithm 4, the value for $\rho$ is chosen such that the
stated function is strongly convex. In addition, in the particular application
of I-BFGS-SC for the minimization of $F$ in (2) that we explore through
numerical experimentation in the following section, we employ the distribution
of terms as described in Section 2.4.

In general, I-BFGS-type methods have no convergence guarantees when employed
to minimize a function that is nonsmooth and nonconvex. However, using the
theoretical guarantees from [42], one can establish a convergence guarantee
for the minimization of an approximate objective if, after some iteration, the
smoothing functions that are employed yield twice-continuously differentiable
functions with Lipschitz Hessians, the smoothing parameters are no longer
modified, the “base point” used for convexification is no longer changed, and
in place of gradients of $f_{i}^{S}(\cdot,\mu_{i})$ one uses gradients of the
convexified terms. In other words, one can establish a convergence guarantee
in a setting in which, after some iteration, the I-BFGS strategy is employed
to minimize a fixed function with strongly convex and sufficiently smooth
terms. The function being minimized in this setting is not the original
nonconvex and nonsmooth objective function, but if the “base point” is close
to a minimizer of the original objective function, then the minimizer of the
employed strongly convex approximation may be close (or closer) to a minimizer
of the original objective.

Under the conditions stated in the previous paragraph (and formalized in the
statement of the theorem below), one can establish the following convergence
guarantee for I-BFGS [42], which we state here for convenience. For simplicity
of notation, we write the guarantee in terms of the minimization of
$\tilde{f}=\tfrac{1}{m}\sum_{i=1}^{m}\tilde{f}_{i}$, where
$\tilde{f}_{i}\approx f_{i}$ for all $i\in\\{1,\dots,m\\}$ such that
$\\{\tilde{f}_{i}\\}_{i=1}^{m}$ satisfies the following assumption. Recall
that we have already established in Lemma 2.5 a set of conditions—including
strong convexity—that ensure the first set of bounds in the assumption. The
final inequality in the assumption can be satisfied if one employs certain
smoothing functions.

######  Assumption 3.1.

There exist $\mu_{\nabla f}\in(0,\infty)$, $L_{\nabla f}\in(0,\infty)$, and
$L_{\nabla^{2}f}\in(0,\infty)$ such that, for all $i\in\\{1,\dots,m\\}$ and
$(\omega,\bar{\omega})$, one has

$\displaystyle\mu_{\nabla
f}\|\omega-\bar{\omega}\|_{2}^{2}\leq(\nabla\tilde{f}_{i}(\omega)-\nabla\tilde{f}_{i}(\bar{\omega}))^{T}(\omega-\bar{\omega})$
$\displaystyle\leq L_{\nabla f}\|\omega-\bar{\omega}\|_{2}^{2}$
$\displaystyle\text{and}\ \
\|\nabla^{2}\tilde{f}_{i}(\omega)-\nabla^{2}\tilde{f}_{i}(\bar{\omega})\|_{2}$
$\displaystyle\leq L_{\nabla^{2}f}\|\omega-\bar{\omega}\|_{2}.$

Since $\tilde{f}$ is strongly convex, it has a unique global minimizer
$\omega_{*}$.

######  Lemma 3.2.

([42, Lemma 1]) Suppose that $i_{k}=(k\ {\rm mod}\ n)+1$ for all
$k\in\mathbb{N}$ and define
$\sigma_{k}:=\max\\{\|z_{k+1,i_{k}}-\omega_{*}\|_{2},\|z_{k,i_{k}}-\omega_{*}\|_{2}\\}$
and $M_{i_{k}}:=\nabla^{2}f_{i_{k}}(\omega_{*})^{-1/2}$ for all
$k\in\mathbb{N}$. If Assumption 3.1 holds and
$\sigma_{k}<m/(3L_{\nabla^{2}f})$, then

$\|B_{k+1,i_{k}}-\nabla^{2}f_{i_{k}}(\omega_{*})\|_{M_{i_{k}}}\leq\left((1-\alpha\theta_{k}^{2})^{1/2}+\alpha_{1}\sigma_{k}\right)\|B_{k,i}-\nabla^{2}f_{i_{k}}(\omega_{*})\|_{M_{i_{k}}}+\alpha_{2}\sigma_{k}$

for some $(\alpha,\alpha_{1},\alpha_{2})\in(0,\infty)^{3}$, where

$\theta_{k}:=\begin{cases}\frac{\|M_{i_{k}}(B_{k,i_{k}}-\nabla^{2}f_{i_{k}}(\omega_{*}))s_{k,i_{k}}\|_{2}}{\|B_{k,i_{k}}-\nabla^{2}f_{i_{k}}(\omega_{*})\|_{M_{i_{k}}}\|M_{i_{k}}^{-1}s_{k,i}\|_{2}}&\text{if
$B_{k,i_{k}}\neq\nabla^{2}f_{i_{k}}(\omega_{*})$}\\\
0&\text{otherwise.}\end{cases}$

######  Lemma 3.3.

([42, Lemma 3]) Suppose that $i_{k}=(k\ {\rm mod}\ n)+1$ for all
$k\in\mathbb{N}$ and define $M_{i}:=\nabla^{2}f_{i}(\omega_{*})^{-1/2}$ for
all $i\in\\{1,\dots,m\\}$. If Assumption 3.1 holds, then, for any $r\in(0,1)$,
there exist $\epsilon(r)\in(0,\infty)$ and $\delta(r)\in(0,\infty)$ such that
having both $\|\omega_{0}-\omega_{*}\|_{2}<\epsilon(r)$ and
$\|B_{0,i}-\nabla^{2}f_{i}(\omega_{*})\|_{M_{i}}<\delta(r)$ for all
$i\in\\{1,\dots,m\\}$ implies that

$\|w_{k}-\omega_{*}\|_{2}\leq
r^{[\frac{k-1}{m}+1]+1}\|\omega_{0}-\omega_{*}\|_{2}.$

In addition, the sequences $\\{\|B_{k,i}\|_{2}\\}$ and
$\\{\|B_{k,i}^{-1}\|_{2}\\}$ are uniformly bounded.

######  Theorem 3.4.

([42, Theorem 7]) Suppose that the conditions in Lemmas 3.2 and 3.3 hold and
define
$\tilde{\omega}_{k}:=argmax_{j\in\\{kn,\dots,kn+n-1\\}}\\{\|\omega_{j}-\omega_{*}\|_{2}\\}$.
Then, $\\{\tilde{\omega}_{k}\\}$ converges to $\omega_{*}$ at a superlinear
rate in the sense that

$\lim_{k\to\infty}\frac{\|\tilde{\omega}_{k+1}-\omega_{*}\|_{2}}{\|\tilde{\omega}_{k}-\omega_{*}\|_{2}}=0.$

In addition, there exists $\\{\zeta_{k}\\}$ such that
$\|w_{k}-\omega_{*}\|_{2}\leq\zeta_{k}$ for all $k\in\mathbb{N}$, where
$\\{\zeta_{k}\\}$ vanishes at a superlinear rate in the sense that

$\lim_{k\to\infty}\frac{\zeta_{k+1}}{\zeta_{k}}=0.$

##  4 Computational Results

We implemented I-BFGS, I-BFGS-DC, I-BFGS-S, I-BFGS-C, and I-BFGS-SC as
described in the previous section to minimize the semi-supervised machine
learning objective function (2) using the data sets that are summarized in
Table 1 below. The data sets that are numbered 1–6 and 9–12 are from [43] and
[14], respectively, while the data sets that are numbered 7–8 were generated
as described in [16]. We ran the experiments using the polyps cluster in the
COR@L Laboratory at Lehigh
University.111https://coral.ise.lehigh.edu/wiki/doku.php/info:coral For the
data sets numbered 1–8, ten-fold cross-validation was performed to split the
data into training and testing sets. For the data sets numbered 9–11, the data
was divided into training and testing sets as described in [14]. For the data
set numbered 12, the data was split into training and testing sets randomly,
as in [3].

Table 1: Data Sets $\\#$ |  Data Set |  $n$ | $m$ |  $\\#$ |  Data Set |  $n$ | $m$
---|---|---|---|---|---|---|---
Training |  Test |  Training |  Test
1 |  Ionosphere |  34 |  315 |  35 |  9 |  a9a |  123 |  32562 |  16561
2 |  Pima |  9 |  690 |  77 |  10 |  w8a |  300 |  49749 |  14951
3 |  Sonar |  60 |  186 |  21 |  11 |  ijcnn1 |  22 |  49990 |  91701
4 |  Diagnostic |  30 |  511 |  57 |  12 |  Covtype |  54 |  464810 |  116202
5 |  Heart |  13 |  243 |  27 |  |  |  |  |
6 |  Cancer |  9 |  614 |  69 |  |  |  |  |
7 |  g50c |  50 |  495 |  55 |  |  |  |  |
8 |  g10n |  10 |  495 |  55 |  |  |  |  |

###  4.1 Parameters of Algorithms

For each data set, the selection of the values $C_{1}$ and $C_{2}$ was
performed as in [3]. In particular, using the training set for each data set,
a training phase was performed that considered the values $C_{1}=10^{i}$ for
each $i\in\\{-1,0,1,2\\}$ and $C_{2}=C_{1}^{-j}$ for each $j\in\\{0,1,2\\}$.
Using these results, the values of $C_{1}$ and $C_{2}$ were set (separately
for each data set) as those yielding the best results for the testing data set
for each data set.

For all variants of I-BFGS, the iteration limit was set as $10^{4}$ for data
sets 1–8 and $2m$ for data sets 9–12. In all cases, the starting point was
generated randomly in the interval $[-5,5]$ for each component of
$(w,b)=(w_{1},...,w_{n},b)$ and the value $c$ was set at $10^{-8}.$ For
I-BFGS-S, I-BFGS-C, and I-BFGS-SC, the parameters $\mu_{0,i}$ for all
$i\in\\{1,\dots,m\\}$, $\sigma$, and $\kappa$ were set as $0.1$, $0.9$, and
$0.5$, respectively. Lastly, for I-BFGS-C and I-BFGS-SC, the $\rho$ values
were set as described Sections 2.3 and 2.4.

For the purposes of comparison, we also implemented the bundle method,
referred to as TSVM-Bundle, with the method-specific parameters as described
in [3, 28]. The iteration limit and starting points were set in the same
manner as for the I-BFGS methods. The subproblems for TSVM-Bundle were solved
using Gurobi [33].

###  4.2 Comparison of Algorithms

Table LABEL:table:results provides the results for all algorithms applied to
each problem instance in terms of average testing error. For each data set, we
consider different percentages of labeled data (indicated in the “% Labeled
Data” column) to show results across a wide spectrum of the balance of terms
in the objective (2). Since ours are only prototype implementations of the
methods, we do not provide detailed time comparisons of the runs. That said,
we mention in passing that the computational time for the I-BFGS methods was
always significant less than the computational times for TSVM-Bundle. For a
glimpse of this fact, we note that the average times for the I-BFGS methods
for data sets 1–8 was around 1.8 seconds and for data sets 9–12 was around 50
seconds, whereas for TSVM-Bundle the average times was on the order of
$10^{3}$ seconds for data sets 1–8 and on the order of $10^{5}$ seconds for
data sets 9–12.

To provide an easier-to-visualize summary of the performance of the
algorithms, we also provide performance profiles [26]. A performance profile
makes use of a performance measure $t_{p,s}$, performance ratio $r_{p,s}$, and
(cumulative) fractional performance $\rho_{s}(\tau)$, where $p$ and $s$ are
indices for problem and algorithm, respectively. For our comparison, we the
relevant improvement of results for the performance measurement (see e.g. [1,
4, 52]), specifically,

$t_{p,s}=1+\frac{f_{p,s}-f^{b}_{p}}{f^{w}_{p}-f^{b}_{p}},$

where $f_{p,s}$ is the objective value of problem $p$ obtained by algorithm
$s$, and $f^{b}_{p}$ and $f^{w}_{p}$ are the best and the worst final
objective values of problem $p$ obtained among all algorithms, respectively.
Generally speaking, a performance profile shows how well each algorithm
performs relative the others. In each profile, a plotted point
$(\tau,\rho_{s}(\tau))$ says that algorithm $s$ finds solutions within $\tau$
times of the best found solution in $\rho_{s}(\tau)*100\%$ of the problems.

The profile in Figure 3 shows that the I-BFGS algorithms all outperform TSVM-
Bundle in general, with I-BFGS-SC yielding the best results overall. Since, as
previously mentioned, the I-BFGS methods require significantly less CPU time,
our results show overall that I-BFGS methods are preferable in this setting.

Table 2: Average testing error of I-BFGS, I-BFGS-DC, I-BFGS-S, I-BFGS-C, I-BFGS-SC and TSVM-Bundle # |  % Labeled Data |  I-BFGS |  I-BFGS-DC |  I-BFGS-S |  I-BFGS-C |  I-BFGS-SC |  TSVM-Bundle
---|---|---|---|---|---|---|---
1 |  10% |  17.71 |  21.14 |  21.14 |  29.43 |  21.14 |  19.14
20% |  16.85 |  14.28 |  22 |  18.29 |  16.86 |  19.71
30% |  14.57 |  4.57 |  13.71 |  14.86 |  15.14 |  13.42
40% |  17.99 |  14 |  16.28 |  15.71 |  18.29 |  16
50% |  13.71 |  13.99 |  12.85 |  12.57 |  12.57 |  14.28
60% |  12.57 |  12.57 |  12 |  12 |  11.71 |  11.42
70% |  11.41 |  13.42 |  13.14 |  13.43 |  12.29 |  13.71
80% |  11.71 |  12 |  12.57 |  12 |  11.43 |  12.57
90% |  12.28 |  11.71 |  10.57 |  11.14 |  11.14 |  11.14
100% |  10.85 |  12.85 |  11.71 |  12.29 |  12.57 |  10.85
2 |  10% |  26.59 |  26.84 |  27.11 |  26.2 |  23.59 |  24.89
20% |  24.75 |  24.23 |  23.71 |  23.97 |  24.23 |  23.7
30% |  22.93 |  22.15 |  23.06 |  23.32 |  22.54 |  22.66
40% |  21.75 |  22.14 |  22.02 |  22.54 |  21.64 |  21.63
50% |  22.02 |  21.49 |  22.02 |  22.81 |  22.01 |  21.75
60% |  20.72 |  20.45 |  21.5 |  21.63 |  21.11 |  21.63
70% |  21.89 |  21.24 |  21.5 |  20.85 |  22.03 |  21.89
80% |  22.15 |  20.85 |  22.68 |  22.29 |  22.02 |  21.63
90% |  22.02 |  20.45 |  22.15 |  22.29 |  21.51 |  21.89
100% |  21.89 |  22.67 |  21.76 |  22.02 |  22.29 |  22.28
3 |  10% |  36.26 |  40.54 |  32.92 |  22.69 |  16.74 |  17.21
20% |  23 |  13.92 |  27.66 |  20.45 |  20.93 |  23.14
30% |  13.16 |  8.76 |  9.11 |  27.4 |  15.93 |  10.16
40% |  19.52 |  8.14 |  17.42 |  13.4 |  15.6 |  29.71
50% |  7.21 |  20.47 |  24.35 |  27.9 |  18.93 |  29.95
60% |  10.07 |  18.47 |  18.78 |  11.5 |  20.31 |  20.59
70% |  26.14 |  21.42 |  13.02 |  13.5 |  15.93 |  17.95
80% |  16.19 |  18 |  20.47 |  15.98 |  26.52 |  24.63
90% |  25.4 |  13.5 |  32.02 |  17.31 |  17.1 |  23.5
100% |  21.26 |  32 |  32 |  31.6 |  29.1 |  37.38
4 |  10% |  14.1 |  7.94 |  14.63 |  11.46 |  12.7 |  11.83
20% |  7.58 |  5.46 |  9.01 |  6.19 |  5.47 |  4.06
30% |  4.48 |  4.4 |  4.41 |  4.59 |  4.76 |  4.06
40% |  4.04 |  3.69 |  4.58 |  3.7 |  2.99 |  2.28
50% |  2.98 |  2.99 |  2.28 |  3.17 |  2.99 |  1.93
60% |  3.58 |  2.98 |  3.33 |  2.99 |  2.28 |  2.28
70% |  2.81 |  3.33 |  2.1 |  2.99 |  3.17 |  2.11
80% |  3.16 |  3.34 |  3.16 |  2.29 |  2.29 |  2.1
90% |  2.1 |  2.45 |  2.63 |  1.76 |  2.46 |  2.45
100% |  2.28 |  2.81 |  2.28 |  3.34 |  2.28 |  2.11
5 |  10% |  28.51 |  20.74 |  23.7 |  24.44 |  25.19 |  25.92
20% |  19.62 |  19.25 |  20.37 |  22.96 |  22.22 |  18.88
30% |  15.92 |  15.18 |  14.81 |  17.78 |  14.44 |  15.18
40% |  16.66 |  13.33 |  16.66 |  15.19 |  15.93 |  16.29
50% |  14.44 |  13.7 |  15.55 |  13.7 |  14.81 |  15.92
60% |  15.55 |  13.7 |  14.44 |  14.81 |  15.19 |  14.44
70% |  14.44 |  14.07 |  13.7 |  14.44 |  14.07 |  14.44
80% |  12.22 |  12.22 |  14.07 |  12.96 |  12.22 |  13.33
90% |  14.07 |  13.33 |  15.18 |  12.59 |  12.96 |  12.96
100% |  14.44 |  14.44 |  14.81 |  14.81 |  15.56 |  14.44
6 |  10% |  3.22 |  3.66 |  2.63 |  3.22 |  2.64 |  2.92
20% |  4.23 |  2.92 |  2.63 |  2.64 |  2.78 |  2.92
30% |  3.36 |  3.21 |  3.21 |  2.63 |  2.78 |  3.5
40% |  3.65 |  3.21 |  3.2 |  2.2 |  2.77 |  3.35
50% |  3.21 |  2.77 |  3.06 |  2.2 |  2.92 |  2.92
60% |  3.35 |  2.77 |  2.63 |  2.34 |  2.92 |  2.77
70% |  3.36 |  3.49 |  2.92 |  2.63 |  2.34 |  2.77
80% |  2.63 |  3.06 |  2.77 |  2.48 |  3.36 |  2.63
90% |  3.21 |  2.92 |  2.48 |  2.05 |  2.63 |  3.07
100% |  2.48 |  2.48 |  3.21 |  2.92 |  3.06 |  2.63
7 |  10% |  30 |  23.63 |  0.54 |  10.91 |  36.91 |  13.45
20% |  16.18 |  29.27 |  28.18 |  19.09 |  14.91 |  28.54
30% |  39.27 |  18 |  9.63 |  38 |  33.45 |  39.45
40% |  33.81 |  27.27 |  19.81 |  30 |  10 |  23.45
50% |  29.81 |  10 |  20.9 |  29.09 |  3.45 |  24
60% |  13.09 |  6.9 |  26 |  12.36 |  24.18 |  30.54
70% |  7.09 |  6.18 |  8.36 |  7.45 |  6.73 |  8.72
80% |  3.63 |  6.72 |  3.63 |  3.64 |  3.82 |  7.27
90% |  4.18 |  4.36 |  4.9 |  5.82 |  4 |  4.72
100% |  4.9 |  1.63 |  6 |  6.18 |  6.73 |  6.72
8 |  10% |  12.18 |  12.54 |  12 |  12.73 |  10.36 |  12.36
20% |  6.9 |  4.72 |  5.27 |  4.91 |  2.73 |  5.09
30% |  4.54 |  3.45 |  4.36 |  3.45 |  2.36 |  6.18
40% |  3.45 |  3.09 |  3.27 |  3.27 |  2.36 |  3.45
50% |  2.36 |  1.45 |  2.54 |  2.18 |  1.27 |  3.09
60% |  1.45 |  1.81 |  1.27 |  1.45 |  0.91 |  1.99
70% |  1.63 |  1.63 |  1.63 |  1.82 |  1.45 |  2
80% |  1.81 |  1.63 |  2.18 |  2.18 |  1.64 |  2.36
90% |  2 |  0.54 |  2 |  1.45 |  1.09 |  3.63
100% |  1.47 |  1.81 |  1.63 |  1.45 |  0.73 |  2.36
9 |  10% |  23.78 |  21.9 |  21.09 |  23.63 |  19.3 |  19.02
20% |  20.89 |  23.63 |  17.23 |  23.54 |  23.6 |  15.86
30% |  20.42 |  34.46 |  18.19 |  17.03 |  23.43 |  24
40% |  17.36 |  20.59 |  20.57 |  23.63 |  17.57 |  17.72
50% |  15.24 |  23.11 |  23.63 |  23.62 |  17.61 |  22.15
60% |  19.93 |  17.49 |  22.28 |  23.62 |  23.63 |  15.77
70% |  20.1 |  21.57 |  18.41 |  21.74 |  23.63 |  20.73
80% |  18.86 |  23.09 |  17.36 |  23.63 |  20.66 |  15.48
90% |  16.3 |  23.58 |  23.63 |  22.12 |  17.37 |  16.32
100% |  17.64 |  23.63 |  23.63 |  23.63 |  23.44 |  15.75
10 |  10% |  8.49 |  48.53 |  12.23 |  7.24 |  6.06 |  9.39
20% |  8.65 |  28.61 |  3.58 |  3.04 |  3.04 |  6.2
30% |  3.12 |  3.04 |  3.04 |  3.04 |  3.04 |  4.32
40% |  2.72 |  3.04 |  3.04 |  3.04 |  3.04 |  3.54
50% |  3.04 |  3.04 |  3.04 |  3.04 |  3.04 |  4.19
60% |  3.04 |  3.04 |  3.04 |  3.04 |  3.04 |  3.63
70% |  3.01 |  3.04 |  3.04 |  3.04 |  3.04 |  3.5
80% |  3.04 |  3.04 |  3.04 |  3.04 |  3.04 |  3.58
90% |  3.04 |  3.04 |  3.04 |  3.04 |  3.04 |  3.34
100% |  3.03 |  3.04 |  3.04 |  3.04 |  3.04 |  3.14
11 |  10% |  9.5 |  9.51 |  9.51 |  9.84 |  9.51 |  9.51
20% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.51
30% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.51
40% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.51
50% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.5
60% |  9.49 |  9.51 |  9.51 |  9.51 |  9.51 |  9
70% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.91
80% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.51
90% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.46
100% |  9.51 |  9.51 |  9.51 |  9.51 |  9.51 |  9.37
12 |  10% |  40.36 |  40.34 |  57.64 |  40.36 |  40.05 |  45.27
20% |  40.36 |  40.36 |  40.36 |  40.36 |  39.63 |  39.56
30% |  39.58 |  58.36 |  40.25 |  40.36 |  39.89 |  54.94
40% |  40.27 |  40.36 |  40.36 |  40.36 |  40.36 |  50.98
50% |  40.34 |  40.36 |  40.36 |  56.89 |  40.36 |  59.24
60% |  40.34 |  40.36 |  40.36 |  57.81 |  40.36 |  55.41
70% |  38.62 |  40.36 |  40.36 |  57.2 |  40.36 |  53.52
80% |  39.87 |  40.36 |  40.36 |  40.36 |  40.36 |  59.65
90% |  40.36 |  40.36 |  51.73 |  40.36 |  33.05 |  59.51
100% |  35.65 |  38.14 |  40.36 |  40.36 |  46.3 |  44.89
Figure 3: Performance profile of algorithms.

##  5 Conclusion

In this paper, we have proposed and demonstrated the numerical performance of
I-BFGS methods for the minimization of a sum of a finite number of functions
that might be nonsmooth and/or nonconvex. In particular, we focused on a
problem that arises in semi-supervised machine learning. Our results show that
variants of I-BFGS that consider a difference-of-convex function formulation,
smoothing, and/or (strongly) convex local approximations outperform a bundle
method that has been designed for this specific problem formulation.

## Disclosure statement

The authors report there are no competing interests to declare.

## Funding

The first author is supported by the Scientific and Technological Research
Council of Türkiye (TUBITAK) with Postdoctoral Research Fellowship Program.

## References

  * [1] M.M. Ali, C. Khompatraporn, and Z.B. Zabinsky, _A numerical evaluation of several stochastic algorithms on selected continuous global optimization test problems_ , Journal of Global Optimization 31 (2005), pp. 635–672.
  * [2] A. Astorino and A. Fuduli, _Nonsmooth optimization techniques for semisupervised classification_ , IEEE Transactions on Pattern Analysis and Machine Intelligence 29 (2007), pp. 2135–2142.
  * [3] A. Astorino and A. Fuduli, _Support vector machine polyhedral separability in semisupervised learning_ , Journal of Optimization Theory and Applications 164 (2015), pp. 1039–1050.
  * [4] V. Beiranvand, W. Hare, and Y. Lucet, _Best practices for comparing optimization algorithms_ , Optimization and Engineering 18 (2017), pp. 815–848.
  * [5] M. Belkin, P. Niyogi, and V. Sindhwani, _Manifold regularization: A geometric framework for learning from labeled and unlabeled examples._ , Journal of Machine Learning Research 7 (2006).
  * [6] K. Bennett and A. Demiriz, _Semi-supervised support vector machines_ , Advances in Neural Information Processing Systems 11 (1998).
  * [7] D.P. Bertsekas, _Incremental gradient, subgradient, and proximal methods for convex optimization: A survey_ , Optimization for Machine Learning 2010 (2011), p. 3.
  * [8] D. Blatt, A.O. Hero, and H. Gauchman, _A convergent incremental gradient method with a constant step size_ , SIAM Journal on Optimization 18 (2007), pp. 29–51.
  * [9] C.G. Broyden, J.E. Dennis, Jr., and J.J. Moré, _On the local and superlinear convergence of quasi-newton methods_ , J. Inst. Math. Applics. 12 (1973), pp. 223–245.
  * [10] C.G. Broyden, _The convergence of a class of double-rank minimization algorithms 1. general considerations_ , IMA Journal of Applied Mathematics 6 (1970), pp. 76–90.
  * [11] R.H. Byrd, G.M. Chin, W. Neveitt, and J. Nocedal, _On the use of stochastic hessian information in optimization methods for machine learning_ , SIAM Journal on Optimization 21 (2011), pp. 977–995.
  * [12] R.H. Byrd, S.L. Hansen, J. Nocedal, and Y. Singer, _A stochastic quasi-newton method for large-scale optimization_ , SIAM Journal on Optimization 26 (2016), pp. 1008–1031.
  * [13] R.H. Byrd, J. Nocedal, and Y.X. Yuan, _Global convergence of a cass of quasi-newton methods on convex problems_ , SIAM Journal on Numerical Analysis 24 (1987), pp. 1171–1190.
  * [14] C.C. Chang and C.J. Lin, _Libsvm: a library for support vector machines_ , ACM Transactions on Intelligent Systems and Technology (TIST) 2 (2011), pp. 1–27.
  * [15] O. Chapelle, V. Sindhwani, and S.S. Keerthi, _Optimization techniques for semi-supervised support vector machines._ , Journal of Machine Learning Research 9 (2008).
  * [16] O. Chapelle and A. Zien, _Semi-supervised classification by low density separation_ , in _International workshop on artificial intelligence and statistics_. PMLR, 2005, pp. 57–64.
  * [17] X. Chen, _Smoothing methods for nonsmooth, nonconvex minimization_ , Mathematical Programming 134 (2012), pp. 71–99.
  * [18] X. Chen, L. Niu, and Y. Yuan, _Optimality conditions and a smoothing trust region Newton method for nonLipschitz optimization_ , SIAM Journal on Optimization 23 (2013), pp. 1528–1552.
  * [19] F.H. Clarke, _Optimization and Nonsmooth Analysis_ , Canadian Mathematical Society Series of Monographs and Advanced Texts, John Wiley & Sons, New York, NY, USA, 1983.
  * [20] R. Collobert, F. Sinz, J. Weston, L. Bottou, and T. Joachims, _Large scale transductive svms._ , Journal of Machine Learning Research 7 (2006).
  * [21] Y. Cui and J.S. Pang, _Modern Nonconvex Nondifferentiable Optimization_ , SIAM, 2021.
  * [22] F.E. Curtis and X. Que, _A quasi-newton algorithm for nonconvex, nonsmooth optimization with global convergence guarantees_ , Mathematical Programming Computation 7 (2015), pp. 399–428.
  * [23] F.E. Curtis, D.P. Robinson, and B. Zhou, _A self-correcting variable-metric algorithm framework for nonsmooth optimization_ , IMA Journal of Numerical Analysis 40 (2020), pp. 1154–1187.
  * [24] D. Davis, D. Drusvyatskiy, K.J. MacPhee, and C. Paquette, _Subgradient methods for sharp weakly convex functions_ , Journal of Optimization Theory and Applications 179 (2018), pp. 962–982.
  * [25] J.E. Dennis and J.J. Moré, _A characterization of superlinear convergence and its application to quasi-newton methods_ , Mathematics of Computation 28 (1974), pp. 549–560.
  * [26] E.D. Dolan and J.J. Moré, _Benchmarking optimization software with performance profiles_ , Mathematical Programming 91 (2002), pp. 201–213.
  * [27] R. Fletcher, _A new approach to variable metric algorithms_ , The Computer Journal 13 (1970), pp. 317–322.
  * [28] A. Fuduli, M. Gaudioso, and G. Giallombardo, _Minimizing nonconvex nonsmooth functions via cutting planes and proximity control_ , SIAM Journal on Optimization 14 (2004), pp. 743–756.
  * [29] D. Goldfarb, _A family of variable-metric methods derived by variational means_ , Mathematics of Computation 24 (1970), pp. 23–26.
  * [30] M. Gurbuzbalaban, A. Ozdaglar, and P.A. Parrilo, _Convergence rate of incremental gradient and incremental newton methods_ , SIAM Journal on Optimization 29 (2019), pp. 2542–2565.
  * [31] M. Gürbüzbalaban, A. Ozdaglar, and P. Parrilo, _A globally convergent incremental newton method_ , Mathematical Programming 151 (2015), pp. 283–313.
  * [32] M. Gurbuzbalaban, A. Ozdaglar, and P.A. Parrilo, _On the convergence rate of incremental aggregated gradient algorithms_ , SIAM Journal on Optimization 27 (2017), pp. 1035–1048.
  * [33] Gurobi Optimization, LLC, _Gurobi Optimizer Reference Manual_ (2022). Available at https://www.gurobi.com.
  * [34] R. Horst and N.V. Thoai, _Dc programming: overview_ , Journal of Optimization Theory and Applications 103 (1999), pp. 1–43.
  * [35] Y. Hu, C.K.W. Yu, and X. Yang, _Incremental quasi-subgradient methods for minimizing the sum of quasi-convex functions_ , Journal of Global Optimization 75 (2019), pp. 1003–1028.
  * [36] H. Iiduka, _Incremental subgradient method for nonsmooth convex optimization with fixed point constraints_ , Optimization Methods and Software 31 (2016), pp. 931–951.
  * [37] T. Joachims, _et al._ , _Transductive inference for text classification using support vector machines_ , in _Icml_ , Vol. 99. 1999, pp. 200–209.
  * [38] K.C. Kiwiel, _Convergence of approximate and incremental subgradient methods for convex optimization_ , SIAM Journal on Optimization 14 (2004), pp. 807–840.
  * [39] A.S. Lewis and M.L. Overton, _Nonsmooth optimization via quasi-newton methods_ , Mathematical Programming 141 (2013), pp. 135–163.
  * [40] D.H. Li and M. Fukushima, _A modified bfgs method and its global convergence in nonconvex minimization_ , Journal of Computational and Applied Mathematics 129 (2001), pp. 15–35.
  * [41] D.H. Li and M. Fukushima, _On the global convergence of the bfgs method for nonconvex unconstrained optimization problems_ , SIAM Journal on Optimization 11 (2001), pp. 1054–1064.
  * [42] A. Mokhtari, M. Eisen, and A. Ribeiro, _Iqn: An incremental quasi-newton method with local superlinear convergence rate_ , SIAM Journal on Optimization 28 (2018), pp. 1670–1698.
  * [43] P.M. Murphy, _Uci repository of machine learning databases_ , ftp:/pub/machine-learning-databaseonics. uci. edu (1994).
  * [44] A. Nedic and D.P. Bertsekas, _Incremental subgradient methods for nondifferentiable optimization_ , SIAM Journal on Optimization 12 (2001), pp. 109–138.
  * [45] J. Nocedal and S. Wright, _Numerical Optimization_ , Springer Science & Business Media, 2006.
  * [46] P. Qi, W. Zhou, and J. Han, _A method for stochastic L-BFGS optimization_ , in _2017 IEEE 2nd International Conference on Cloud Computing and Big Data Analysis (ICCCBDA)_. IEEE, 2017, pp. 156–160.
  * [47] S.S. Ram, A. Nedić, and V.V. Veeravalli, _Incremental stochastic subgradient algorithms for convex optimization_ , SIAM Journal on Optimization 20 (2009), pp. 691–717.
  * [48] D.F. Shanno, _Conditioning of quasi-newton methods for function minimization_ , Mathematics of Computation 24 (1970), pp. 647–656.
  * [49] V. Sindhwani and S.S. Keerthi, _Large scale semi-supervised linear SVMs_ , in _Proceedings of the 29th annual international ACM SIGIR conference on Research and development in information retrieval_. 2006, pp. 477–484.
  * [50] M.V. Solodov, _Incremental gradient algorithms with stepsizes bounded away from zero_ , Computational Optimization and Applications 11 (1998), pp. 23–35.
  * [51] J.E. Van Engelen and H.H. Hoos, _A survey on semi-supervised learning_ , Machine Learning 109 (2020), pp. 373–440.
  * [52] A.I.F. Vaz and L.N. Vicente, _A particle swarm pattern search method for bound constrained global optimization_ , Journal of Global Optimization 39 (2007), pp. 197–219.
  * [53] J.P. Vial, _Strong and weak convexity of sets and functions_ , Mathematics of Operations Research 8 (1983), pp. 231–259.
  * [54] H.T. Wai, W. Shi, C.A. Uribe, A. Nedić, and A. Scaglione, _Accelerating incremental gradient optimization with curvature information_ , Computational Optimization and Applications 76 (2020), pp. 347–380.
  * [55] R. Zhao, W.B. Haskell, and V.Y. Tan, _Stochastic l-bfgs: Improved convergence rates and practical acceleration strategies_ , IEEE Transactions on Signal Processing 66 (2017), pp. 1155–1169.
  * [56] X. Zhou and M. Belkin, _Semi-supervised learning_ , in _Academic press Library in signal processing_ , Vol. 1, Elsevier, 2014, pp. 1239–1269.