\section{Introduction}
\label{sec:Introduction}
Statistical analysis using (generalized) self-concordant (SC) functions as a loss function is gaining increasing attention in the machine learning community \cite{Bac10,Owe13,OstBac18,MarBacOstRu19}.  This class of loss functions allows to obtain faster statistical rates akin to least-squares. At the same time, the minimization of empirical risk in this setting is a challenging optimization problem in high dimensions. The reason is that for large-scale problems usually a \ac{FOM} is the method of choice. Yet, self-concordant functions are usually not strongly convex and do not have Lipschitz continuous gradients. Hence, unless we can verify a relative smoothness condition \cite{BauBolTeb17,LuFreNes18},  the cornerstone assumptions for \ac{FOM} do not apply when minimizing \ac{SC} functions. Further, many machine learning problems involve some sparsity constraint usually admitting a \ac{LMO} \cite{Jag13}, meaning that for a reasonable computational price one can minimize a linear function over the feasible set; Examples include $\ell_{1}$ regularization, Spectrahedron constraints \cite{Haz08}, and ordered weighted $\ell_1$ regularization \cite{FigZen14}. In such settings an attractive method of choice is the Frank-Wolfe (\ac{FW}) method \cite{FraWol56}, also known as the \ac{CG}. The modern convergence analysis of \ac{FW}, starting with \cite{Jag13}, relies on bounded curvature of the objective function $f$. A sufficient condition for this is that the objective function has a Lipschitz continuous gradient over the domain of interest. In general, \ac{SC} functions have unbounded curvature and are not strongly convex. In fact, to the best of our knowledge, there is no existing convergence guarantee for \ac{FW} when minimizing a general \ac{SC} function. Given the plethora of examples of \ac{SC} loss functions, and the importance of \ac{FW} in machine learning, we focus on the question of convergence and complexity of \ac{FW} for minimizing a SC function $f$ over a set $\mathcal{X}$ with \ac{LMO}. 

\textbf{Our Contributions. } In this paper we demonstrate that \ac{FW} indeed works when minimizing a \ac{SC} function over a compact convex domain. However, since the standard arguments for \ac{FW} are not applicable, we have to develop some new ideas. Section \ref{sec:sublinear} constructs new adaptive variants of the \ac{FW} method, whose main innovation is the construction of new step-size policies ensuring global convergence and $\mathcal{O}(1/k)$ convergence rates. One of the keys to our results is to develop step-sizes which ensure monotonicity of the method, leading to Lipschitz smoothness and strong convexity of the objective on a level set. While the Lipschitz smoothness on a level set is used in the analysis, the algorithms do not require knowledge of the Lipschitz constant but rather rely only on basic properties of \ac{SC} functions.
Section \ref{sec:fast} gives a \emph{linear} convergence result for \ac{FW}, under the assumption that we have access to the strong convexity parameter on the level set and a \ac{LLOO}, first constructed in \citet{GarHaz16}. 

\textbf{Previous work. } Despite its great scalability properties, \ac{FW} is plagued by slow convergence rates. \citet{FraWol56} showed an $\mathcal{O}(1/k)$ convergence rate of the function values to the optimal value for quadratic programming problems. Later on \citet{LevPol66} generalized this result to general functions and arbitrary compact convex sets, and \citet{CanCul68} proved that this rate is actually tight (see also \citet{Lan13}). However, a standard result in convex optimization states that projected gradient descent converges linearly to the minimizer of a strongly convex and Lipschitz (i.e. well-conditioned) smooth function. Subsequently, an important line of research emerged with the aim to accelerate \ac{FW}. \citet{GueMar86} obtained linear convergence rates in well conditioned problems over polytopes under the a-priori assumption that the solution lies in the relative interior of the feasible set, and the rate of convergence explicitly depends on the distance the solution is away from the boundary (see also \citet{EpeFreu00,BecTeb04}). If no a-priori information on the location of the solution is available, there are essentially two known twists of the vanilla \ac{FW} to boost the convergence rates. One twist is to modify the search directions via \emph{corrective} or \emph{away} search directions \cite{JagLac15,FreGriMaz17,BecSht17,GutPen20,PenRod18}. These methods require however more fine grained oracles which are not in general available for \ac{SC} functions\footnote{The away step version of \ac{FW} needs a vertex oracle. Vertices need not be in the domain of a \ac{SC} function (i.e. the origin in Example \ref{example}).}. The alternative twist is to restrict the \ac{LMO} in order to gain a more powerful local approximating model of the objective function. This last strategy has been used in several recent papers \cite{HarJudNem15,GarHaz16}. A unified analysis for different settings was recently proposed in \cite{PedNegAskJag20}, which inspired one of the variants of our methods.

None of these references explicitly studied the important case of \ac{SC} minimization.
\citet{Nes18CG} requires the gradient of the objective function to be H\"older continuous. Implicitly it is assumed that the feasible set $\mathcal{X}$ is contained in the domain of the objective function $f$, an assumption we do not make, and is also not satisfied in important applications (e.g. $0\in\mathcal{X}$ in the Poisson inverse problem below, but $0\notin\dom f$). Specialized to solving a quadratic Poisson inverse problem in phase retrieval, \citet{Odor16} provided a globally convergent \ac{FW} method. They analyzed this problem via techniques developed in \cite{Nes18CG} and with step-size specific to the application. We instead give a unified general analysis valid for \emph{all} \ac{SC} functions. Given the importance of \ac{SC} minimization in machine learning, a large volume of work around proximal Newton and proximal gradient methods have been developed in a series of papers by \cite{CevKyrTra13,CevKyrTra14,CevKyrTra15}. This paper complements this line of research by providing a \ac{FW} analysis of \ac{SC} optimization problems. Recently, \cite{LiuCevTra20} coupled a \ac{LMO} with a second-order model of the \ac{SC} objective function and obtained similar results to ours. Similar to \cite{LiuCevTra20}, the recent paper \cite{carderera2020second-order} studied a Newton-\ac{FW} approach, but for well-conditioned functions. 

\section{Preliminaries}
\label{sec:prelims}
We consider the minimization problem
\begin{equation}\label{eq:P
\min_{x\in\mathcal{X}}f(x),
\end{equation}
where $\mathcal{X}$ is a compact convex set in $\mathbb{R}^n$ with nonempty interior. We assume that $\dom f\cap\mathcal{X}\neq\emptyset$ and the set $\mathcal{X}$ is represented by an \ac{LMO} returning at a point $x$ the target vector 
\begin{equation}\label{eq:s}
s(x)\triangleq\arg\min_{s\in\mathcal{X}}\inner{\nabla f(x),s}.
\end{equation} 
In case of multiple solutions of \eqref{eq:s} we assume that ties are broken by some arbitrary mechanism. The function $f$ is assumed to be self-concordant (\ac{SC}) \cite{NesIPM94} in the following sense: 
Let $f\in{\mathbf{C}}^{3}(\dom f)$ be a closed convex function with open domain $\dom f\triangleq\{x\in\mathbb{R}^n\vert f(x)<\infty\}$.  For a fixed $x\in\dom f$ and direction $u\in\mathbb{R}^n$, define $\phi(x;t)\triangleq f(x+tu)$ for $t\in\dom\phi(x;\cdot)$.
\begin{definition}\label{def:SC}
A proper closed convex function $f:\mathbb{R}^n\to(-\infty,\infty]$ with open domain $\dom f$ is self-concordant (SC) with parameter $M>0$ (i.e. $f\in\mathcal{F}_{M}$) iff 
\begin{equation}
\abs{\phi'''(x,0)}\leq M\phi''(x,0)^{3/2}.
\end{equation}
\end{definition}
\ac{SC} functions were originally introduced in the context of interior-point methods \cite{NesIPM94}, but recently have received significant attention in machine learning as many problems in supervised learning and empirical risk minimization involve loss functions which are \ac{SC}. We give a (non-exhaustive) list of examples next. Note that, unlike classical settings, $f$ in \eqref{eq:P} is not assumed to be a self-concordant barrier for $\mathcal{X}$.

\subsection{Examples}
\label{sec:examples}

\textbf{Poisson Inverse Problem. }
Consider the low-light imaging problem in signal processing, where the imaging data is collected by counting photons hitting a detector over time. In this setting, we wish to accurately reconstruct an image in low-light, which leads to noisy measurements due to low photon count. Assume that the data-generating process of the observations follows as Poisson distribution 
\[
p(y\vert Wx^{\ast})=\prod_{i=1}^{m}\frac{(w_{i}^{\top}x^{\ast})^{y_{i}}}{y_{i}!}\exp(-w_{i}^{\top}x^{\ast}),
\]
where $x^{\ast}\in\mathbb{R}^n$ is the true image, $W$ is a linear operator that projects the scene onto observations, $w_{i}$ is the $i$-th row of $W$ and $y\in\mathbb{N}^{m}$ is the vector of observed photon counts. The maximum likelihood formulation of the Poisson inverse problem \cite{HarMarWil11} under sparsity constraints leads to the optimization problem
\begin{equation}\label{eq:Poisson}
\min_{x\in\mathbb{R}^n_{+};\norm{x}_{1}\leq M} \{f(x)=\sum_{i=1}^{m}w_{i}^{\top}x-\sum_{i=1}^{m}y_{i}\ln(w_{i}^{\top}x)\}.
\end{equation}
Setting $\varphi_{i}(t)\triangleq t-y_{i}\ln(t)$ for $t>0$, we observe that $f(x)= \sum_{i=1}^{m}\varphi_{i}(w_{i}^{\top}\theta)$. Since each individual function $\varphi_{i}$ is \ac{SC} with parameter $M_{\varphi_{i}}\triangleq \frac{2}{\sqrt{y_{i}}}$, general rules of \ac{SC} calculus shows that $f$ is \ac{SC} with domain $\bigcap_{i=1}^{m}\{x\in\mathbb{R}^n\vert w_{i}^{\top}x>0\},$ and parameter $M\triangleq\max_{1\leq i\leq m}\frac{2}{\sqrt{y_{i}}}$ (\citet{Nes18}, Thm. 5.1.1).

\textbf{Learning Gaussian Markov random fields. }
We consider learning a Gaussian graphical random field of $p$ nodes/variables from a data set $\{\phi_{1},\ldots,\phi_{N}\}$ \cite{SpeKli86,CevKyrTra13}. Each random vector $\phi_{j}$ is an i.i.d realization from a $p$-dimensional Gaussian distribution with mean $\mu$ and covariance matrix $\Sigma$. Let $\Theta=\Sigma^{-1}$ be the precision matrix. To satisfy conditional dependencies between the random variables, $\Theta$ must have zero in $\Theta_{ij}$ if $i$ and $j$ are not connected in the underlying dependency graph. To learn the graphical model subject to sparsity constraints, we minimize the loss function 
$$
f(x)\triangleq -\ln\det(\mat(x))+\tr(\hat{\Sigma}\mat(x))
$$
over the $\ell_{1}$-ball $\mathcal{X}\triangleq \{x\in\mathbb{R}^n\vert \mat(x)\succ 0,\norm{x}_{1}\leq R\}$, where $n=p^{2}$ and $\mat(x)\in\mathbb{R}^{p\times p}_{\text{sym}}$ is the $p\times p$ symmetric matrix constructed from the $n=p^{2}$-dimensional vector $x$. $f$ is \ac{SC} with $M=2$.

\textbf{Logistic Regression. }
Consider the regularized logistic regression problem 
\begin{equation}\label{eq:Logistic}
f(x)=\frac{1}{N}\sum_{i=1}^{N}\ell(y_{i}(\inner{\phi_{i},x}+\mu))+\frac{\gamma}{2}\norm{x}^{2}_{2}
\end{equation}
where $\ell(t)\triangleq \log(1+e^{-t})$ is the logistic loss, $\mu$ is a given intercept, $y_{i}\in\{-1,1\}$ and $\phi_{i}\in\mathbb{R}^n$ are given as input data for $i=1,2,\ldots,N$, and $\gamma>0$ is a given regularization parameter. According to Prop. 5 in \citet{SunTran18}, this functions is \ac{SC} with parameter $M=\frac{1}{\sqrt{\gamma}}\max\{\norm{\phi_{i}}_{2}\vert 1\leq i\leq n\}$. To promote sparsity, we minimize \eqref{eq:Logistic} over the $\ell_{1}$-ball $\mathcal{X}=\{x\in\mathbb{R}^n\vert \norm{x}_{1}\leq R\}$.  The resulting \ac{SC} optimization problem \eqref{eq:P} is the constrained formulation of the elastic-net  regularization. 

\textbf{Portfolio Optimization. } In this problem there are $n$ assets with returns $r_{t}\in\mathbb{R}^{n}_{+}$ in period $t$ of the investment horizon. The goal is to minimize the utility function $f(x)$ of the investor by choosing the weights of the assets in the portfolio. Our task is to design a portfolio $x$ solving the problem 
\begin{equation}\label{eq:Portfolio}
\min_x \left\{f(x)=-\sum_{t=1}^{T}\ln(r_{t}^{\top}x) : x_{i}\geq 0,\sum_{i=1}^{n}x_{i}=1 \right\}.
\end{equation}
Note that this problem can be cast into an online optimization model \cite{AroHaz06}. It can be seen that $f\in\mathcal{F}_{2}$. We remark that this problem appears also in the classical universal prediction problem in information theory and online learning \cite{MerFed98}.

\subsection{The Limits of Standard \ac{FW}}
In this section we explain why the standard analysis of \ac{FW} does not work for \ac{SC} functions. For this, we need to briefly recapitulate the main steps in the modern analysis initiated in \citet{Cla10}, and later on refined by \citet{Haz08} for matrix problems over a spactrahedron, and finally concluded in \citet{Jag13} for general convex compact domains.

The typical merit function employed in \ac{FW} algorithms is the dual gap function 
\begin{equation}\label{eq:gap}
\gap(x)\triangleq\max_{s\in\mathcal{X}}\inner{\nabla f(x),x-s}.
\end{equation}
It is easy to see that $\gap(x)\geq 0$ for all $x\in\mathcal{X}$, with equality if and only if $x$ is a solution to \eqref{eq:P}.
\begin{algorithm}[t]
 \caption{Standard Frank-Wolfe method}
 \label{alg:FW}
 \begin{algorithmic}
\STATE {\bfseries Input: } $x^{0}\in\dom f\cap \mathcal{X}$ initial state, Step size policy $(\alpha_{k})_{k\geq 0}$ (either $\alpha_{k}=\frac{2}{k+2}$, or via line-search); $\varepsilon>0$ tolerance level
\FOR{$k=0,1,\ldots$}
 \IF{$\gap(x^{k})>\varepsilon$}
  \STATE Obtain $s^{k}=s(x^{k})$
  \STATE Set $x^{k+1}=x^{k}+\alpha_{k}(s^{k}-x^{k})$
\ENDIF
\ENDFOR
\end{algorithmic}
\end{algorithm}
The convergence analysis of Algorithm \ref{alg:FW} relies on the \emph{curvature constant} \cite{Jag13}
\begin{align*}
C_{f}\triangleq \sup_{x,s\in\mathcal{X},\gamma\in[0,1]} \frac{2}{\gamma^{2}} D_{f}(x+\gamma(s-x),x),
\end{align*}
where we define the \emph{Bregman divergence} of the smooth convex function $f$ as 
\begin{align*}
D_{f}(y,x)\triangleq f(y)-f(x)-\inner{\nabla f(x),y-x}
\end{align*}
for all $x,y\in\dom f$. Under the assumption that $C_{f}<\infty$, \citet{Jag13} proved sublinear $O(1/k)$ rate of convergence in terms of the dual gap function $\gap(\cdot)$. Unfortunately, minimizing a \ac{SC} function over a compact set does not necessarily give us a finite curvature constant, as the following example illustrates.

\begin{example}\label{example}
Consider the function $f(x,y)=-\ln(x)-\ln(y)$ considered over the set $\mathcal{X}=\{(x,y)\in\mathbb{R}^{2}_{+}\vert x+y=1\}$. This function is the standard \ac{SC} barrier for the positive orthant (the log-barrier) and its Bregman divergence is easily calculated as 
\begin{align*}
D_{f}(u,v)=\sum_{i=1}^{2}\left[-\ln\left(\frac{u_{i}}{v_{i}}\right)+\frac{u_{i}}{v_{i}}-1\right]
\end{align*}
where $u=(u_{1},u_{2}),v=(v_{1},v_{2})$. We see that $C_{f}=\infty$. Neither the function $f$, nor its gradient, is Lipschitz continuous over the set of interest. Moreover, if we start from $u^0=(1/4,3/4)$, using the standard step-size $\alpha_k=2/(k+2)$, we have $\alpha_0=1$ and $u^1=s(u^0)=(1,0) \notin \dom f$ and the method fails.
\end{example}
We point out that Example \ref{example} is representative for the class of optimization problems of interest in this work since the logarithm is the prime example for a self-concordant function. It is thus clear that the standard analysis based on finite curvature estimates (or, as a particular case, Lipschitz continuous gradient) cannot be applied to analyze \ac{FW} when applied to \eqref{eq:P}.

\subsection{Basic Estimates}
From \citet{Nes18}, Thm. 5.1.6, we know that if $\dom f$ contains no lines, then $\nabla^{2}f(x)\succ 0$ for all $x\in\dom f$ ($f$ is \emph{non-degenerate}). We shall assume that $f$ is non-degenerate. Define the local norm of $u\in\mathbb{R}^n$ at $x\in\dom f$ as $\norm{u}_{x}\triangleq \sqrt{\inner{\nabla^{2}f(x)u,u}}$. The dual norm at $x\in\dom f$ is defined as $\norm{u}^{\ast}_{x}\triangleq\sqrt{\inner{[\nabla^{2}f(x)]^{-1}u,u}}$. Given $f\in\mathcal{F}_{M}$, we define the distance-like function 
$$
\mathsf{d}(x,y)\triangleq\frac{M}{2}\norm{y-x}_{x},
$$
as well as 
$$
\omega(t)\triangleq t-\ln(1+t),\text{ and }\omega_{\ast}(t)\triangleq-t-\ln(1-t).
$$
It is not hard to verify that $\omega(t)\geq 0$ for all $t>-1$ and $\omega_{\ast}(t)\geq 0$ for every $t<1$; Moreover, $\omega$ and $\omega_{\ast}$ are increasing and strictly convex in $[0,\infty)$ and $[0,1)$, respectively. For all $x,\tilde{x}\in\dom f$ Thms. 5.1.8 and 5.1.9 in \citet{Nes18} state that
\begin{align}
\label{eq:down}
f(\tilde{x})&\geq f(x)+\inner{\nabla f(x),\tilde{x}-x}+\frac{4}{M^{2}}\omega\left(\mathsf{d}(x,\tilde{x})\right)\\
f(\tilde{x})&\leq f(x)+\inner{\nabla f(x),\tilde{x}-x}+\frac{4}{M^{2}}\omega_{\ast}\left(\mathsf{d}(x,\tilde{x})\right)
\label{eq:up}
\end{align}
where in the latter inequality we assume $\mathsf{d}(x,\tilde{x})<1$.\\
An attractive feature of \ac{SC} functions, in particular from the point of view of \ac{FW} algorithms \cite{Jag13}, is that self-concordance is invariant under affine transformations and re-parametrizations of the domain. See Section A of the supplementary materials for a precise statement of this fact, and further properties of \ac{SC} functions.
Affine invariance allows us to cover a very broad class of composite convex optimization problems of the form
\begin{equation}\label{eq:ERM}
\min_{x\in\mathcal{X}}\{f(x)\triangleq g(\mathbf{E} x)+\inner{q,x}\},
\end{equation}
where $g$ is \ac{SC}, and $\mathbf{E}\in\mathbb{R}^{m\times n},q\in\mathbb{R}^n$. This formulation covers empirical risk-minimization with a convex norm-like regularization, as explained in Section \ref{sec:examples}, among many others.

We have the following existence result for solutions of problem \eqref{eq:P}. For a proof see Section A in the supplementary materials.
\begin{proposition}\label{prop:existence}
Suppose there exists $x\in\dom f\cap\mathcal{X}$ such that $\norm{\nabla f(x)}_{x}^{\ast}\leq \frac{2}{M}$. Then \eqref{eq:P} admits a unique solution.
\end{proposition}
We note that there can be at most one solution to problem \eqref{eq:P} because $\nabla^{2}f\succ 0$ on $\dom f\cap\mathcal{X}$. Henceforth, we will assume that \eqref{eq:P} admits a solution $x^{\ast}$. 
 
\section{\ac{FW} works for Self-Concordant Functions}
\label{sec:sublinear}
This section gives a high-level overview on the necessary modifications of \ac{FW} in order to derive provably convergent schemes for minimizing an \ac{SC} function over a compact convex set $\mathcal{X}$. We propose two new step-size policies, both derived from different approximating local models guaranteeing a version of the celebrated descent lemma. To explain our approach, let us define for all $x\in\dom f$
\begin{equation}\label{eq:e}
\mathtt{e}(x)\triangleq\mathsf{d}(x,s(x))=\frac{M}{2}\norm{s(x)-x}_{x}.
\end{equation}
The first step-size rule we analyze, called Variant 1 (V1), is based on the local sufficient decrease property \eqref{eq:up}. In particular, the step-size policy in V1 minimizes at each iteration a local model in the r.h.s. of \eqref{eq:up}, and reads explicitly as 
\begin{equation}
\alpha_{k}^{V1}\triangleq \min\left\{1,\frac{\gap(x^{k})}{\mathtt{e}(x^{k})(\gap(x^{k})+\frac{4}{M^{2}}\mathtt{e}(x^{k}))}\right\}.
\end{equation}
Observe, that this step-size rule needs to compute the gap function and the local distance $\mathtt{e}(x^{k})$. For evaluating this local distance, it is not necessary to evaluate the full Hessian matrix $\nabla^{2}f(x)$, but rather only the matrix-vector product $\nabla^{2}f(x)(s^{k}-x^{k})$. In many situations of interests, this vector can be computed easily. A particular appealing instance where this applies is the generalized linear model $f(x)=\sum_{i=1}^{N}f_{i}(\inner{a_i,x})$, each function $f_{i}$ being \ac{SC}. In this case, the Hessian matrix is a sum of rank-one matrices, Hessian-vector product amounts to vector-vector and vector-scalar products, and the cost of calculating $\mathtt{e}(x^{k})$ is of the same order as the cost of calculating the gradient.

Variant 2 (V2) is based on a backtracking procedure in the spirit of \citet{PedNegAskJag20}. Specifically, V2 employs the quadratic model 
\begin{equation}\label{eq:Q}
Q(x^{k},t,\mu)\triangleq f(x^{k})-t\gap(x^{k})+\frac{t^{2}\mu}{2}\norm{s(x^{k})-x^{k}}^{2}_{2}
\end{equation}
If $\mu \geq \mathcal{L}_{k}$, where $\mathcal{L}_{k}$ is a local gradient Lipschitz constant estimate, $Q(x^{k},t,\mu)$ is an upper approximation for the objective. This approximation has the advantage that its minimum over $t\in[0,1]$ can be computed in closed form, yielding the step size
\begin{equation}
\alpha_{k}^{V2}\triangleq \min\left\{1,\frac{\gap(x^{k})}{\mathcal{L}_{k}\norm{s(x^{k})-x^{k}}_{2}^{2}}\right\}.
\end{equation}
The quality of this quadratic approximation will depend on the local Lipschitz estimate $\mathcal{L}_{k}$. This parameter needs to be carefully selected to ensure convergence while keeping the number of function evaluations to a minimum. This is achieved through Algorithm \ref{alg:backtrack}, which mimics the backtracking procedure of \citet{PedNegAskJag20}.

\begin{algorithm}
 \caption{Adaptive Frank-Wolfe method for \ac{SC} functions}
 \label{alg:AFW} 
 \begin{algorithmic}
\STATE {\bfseries Input: } $x^{0}\in\dom f\cap \mathcal{X}$ initial state, $\beta\in(0,1)$, Tolerance $\varepsilon>0$.
\FOR{$k=1,\ldots$}
 \IF{$\gap(x^{k})>\varepsilon$}
  \STATE Obtain $s^{k}=s(x^{k})$ and set $v^{k}=s^{k}-x^{k}$.
  \STATE Obtain $\mathtt{e}^{k}=\frac{M}{2}\norm{v^{k}}_{x^{k}}$
  \IF{Variant 1:}
  \STATE Set $\mathtt{t}_{k}=\frac{\gap(x^{k})}{\mathtt{e}^{k}(\gap(x^{k})+\frac{4}{M^{2}}\mathtt{e}^{k})}$
  \STATE Set $\alpha_{k}=\min\{1,\mathtt{t}_{k}\}$
  \ENDIF
  \IF{Variant 2: }
  \STATE Set $(\alpha_{k},\mathcal{L}_{k})=\mathtt{step}(f,v^{k},x^{k},\gap(x^{k}),\mathcal{L}_{k-1})$ 
  \ENDIF
  \STATE Set $x^{k+1}=x^{k}+\alpha_{k}v^{k}$
\ENDIF
\ENDFOR
\end{algorithmic}
\end{algorithm}
The backtracking performed in this subroutine defines a candidate step size $\alpha$ and checks whether the sufficient decrease condition
\begin{equation}\label{eq:SD}
f(x^{k}+\alpha(s^{k}-x^{k}))<Q(x^{k},\alpha,\mu)
\end{equation}
is satisfied. If not, then we increase the proposed Lipschitz estimate to $\gamma_{u}\mathcal{L},\gamma_{u}>1,$ and repeat. As will be shown in the analysis of this method, this described backtracking procedure is guaranteed to stop in finite steps. 
\begin{algorithm}
 \caption{Function $\mathtt{step}(f,v,x,g,\mathcal{L})$}
 \label{alg:backtrack}
\begin{algorithmic}
    \STATE Choose $\gamma_{u}>1,\gamma_{d}<1$
    \STATE Choose $\mu\in[\gamma_{d}\mathcal{L},\mathcal{L}]$
    \STATE $\alpha=\min\{\frac{g}{\mathcal{L}\norm{d}^{2}_{2}},1\}$
    \IF{$f(x+\alpha v)>Q(x,\alpha,\mu)$}
   \STATE $\mu=\gamma_{u}\mu$
    \STATE $\alpha=\min\{\frac{g}{\mu\norm{v}^{2}_{2}},1\}$
    \ENDIF
    \STATE Return $\alpha,\mu$
\end{algorithmic}
\end{algorithm}
\subsection{Analysis of Variant 1}
Section B in the supplementary materials contains complete proofs of all the derivations and claims made in this section.

For $x\in\dom f$, define the target vector $s(x)$ as in \eqref{eq:s}, and $\gap(x)$ as in \eqref{eq:gap}. Given $x\in\mathcal{X}$ and $t>0$, set $x^{+}_{t}\triangleq x+t(s(x)-x)$. Assume that $\mathtt{e}(x)\neq 0$. By construction, 
 \begin{align*}
\mathsf{d}(x,x^{+}_{t})=\frac{tM}{2}\norm{s(x)-x}_{x}=t\mathtt{e}(x)<1,
 \end{align*}
 iff $t<1/\mathtt{e}(x)$. Choosing $t\in(0,1/\mathtt{e}(x))$, we conclude from \eqref{eq:up} 
 \begin{align*}
 f(x^{+}_{t})&\leq f(x)+\inner{\nabla f(x),x^{+}_{t}-x}+\frac{4}{M^{2}}\omega_{\ast}(t\mathtt{e}(x))\\
 &= f(x)-t\gap(x)+\frac{4}{M^{2}}\omega_{\ast}(t\mathtt{e}(x))
 \end{align*}
This shows that when minimizing an \ac{SC} function, we can search for a step size $\alpha_{k}$ which minimizes the r.h.s. of the previous inequality by maximizing
 \begin{equation}\label{eq:eta}
 \eta_{x}(t)\triangleq t\gap(x)-\frac{4}{M^{2}}\omega_{\ast}(t\mathtt{e}(x)),
 \end{equation} 
for  $t\in(0,1/\mathtt{e}(x)).$ As shown in Section B.2 of the supplementary materials, this function attains a unique maximum at the value 
\begin{equation}\label{eq:t}
 \mathtt{t}(x)\triangleq \frac{\gap(x)}{\mathtt{e}(x)(\gap(x)+\frac{4}{M^{2}}\mathtt{e}(x))}\equiv\frac{\gamma(x)}{\mathtt{e}(x)}.
 \end{equation}

We then construct the step size sequence $(\alpha_{k})_{k\geq 0}$ by setting $\alpha_{k}=\min\{1,\mathtt{t}(x^{k})\}$. The convexity of $\mathcal{X}$ and the fact that $\alpha_{k}\mathtt{e}(x^{k})<1$ guarantee that $(x^{k})_{k\geq 0}\subset\dom f\cap\mathcal{X}$. Thus, at each iteration, we reduce the objective function value by at least the quantity $\Delta_{k}\equiv\eta_{x^{k}}(\alpha_{k})$, so that $f(x^{k+1})\leq f(x^{k})-\Delta_{k}<f(x^{k}).$ Hence, we see that the method induces a sequence of function values $(f(x^{k}))_{k\geq 0}$ which is monotonically decreasing by at least the amount $\Delta_{k}$ in each iteration. Equivalently,
\[
 (x^{k})_{k\geq 0}\subset\mathcal{S}(x^{0})\triangleq\{x\in\dom f\cap\mathcal{X}\vert f(x)\leq f(x^{0})\}.
\] 
Furthermore, $\sum_{k\geq 0}\Delta_{k}<\infty$, and hence the sequence $\left(\Delta_{k}\right)_{k\geq 0}$ converges to 0, and 
 \begin{align*}
 \min_{0\leq k<K}\Delta_{k}\leq\frac{1}{K}(f(x^{0})-f^{\ast}).
 \end{align*}

It follows from Lemma B.1 in the supplementary materials that $\frac{4}{M^{2}}\omega(\mathsf{d}(x^{\ast},x))\leq f(x)-f^{\ast}$ for all $x \in \dom f$. Consequently, 
 \begin{align*}
 \mathcal{S}(x^{0})\subseteq\left\{x\in\dom f\vert\omega(\mathsf{d}(x^{\ast},x))\leq \frac{M^{2}}{4}(f(x^{0})-f^{\ast})\right\}.
 \end{align*}
Since $\omega(\cdot)$ is continuous and increasing, $\mathcal{S}(x^{0})$ is a closed, bounded and convex set. Accordingly, defining by $\lambda_{\max}(\nabla^{2}f(x))$ and $\lambda_{\min}(\nabla^{2}f(x))$ the largest, respectively smallest, eigenvalue of the Hessian $\nabla^{2}f(x)$, the numbers 
\begin{align*}
&L_{\nabla f}\triangleq\max_{x\in\mathcal{S}(x^{0})}\lambda_{\max}(\nabla^{2}f(x)),\text{ and } \\
&\sigma_{f}\triangleq \min_{x\in\mathcal{S}(x^{0})}\lambda_{\min}(\nabla^{2}f(x)),
 \end{align*}
 are well defined and finite. In particular, for all $x\in\mathcal{S}(x^{0})$ we have a restricted strong convexity of the function $f$, in the sense that 
 \begin{equation}\label{eq:gbound}
f(x)-f(x^{\ast})\geq \frac{\sigma_{f}}{6}\norm{x-x^{\ast}}_{2}^{2}.
 \end{equation}



 In terms of these quantities, we can bound the sequence $(\mathtt{e}^{k})_{k\geq 0}$, defined as $\mathtt{e}^{k}\triangleq \mathtt{e}(x^{k})$, as
 \begin{equation}\label{eq:normbound}
\frac{M\sqrt{\sigma_{f}}}{2}\norm{s^{k}-x^{k}}_{2}\leq \mathtt{e}^{k}\leq \frac{M\sqrt{L_{\nabla f}}}{2}\norm{s^{k}-x^{k}}_{2}.
 \end{equation}


In order to derive convergence rates, we need to lower bound the per-iteration decrease in the objective function. A detailed analysis, given in Section B of the supplementary materials, reveals that we can construct a useful minorizing sequence for the per-iteration function decrease $(\Delta_{k})_{k\geq 0}$ as 
\begin{equation}\label{eq:lowDelta}
\Delta_{k}\geq \min\{\mathtt{a}\gap(x^{k}),\mathtt{b}\gap(x^{k})^{2}\},
\end{equation}
where $\mathtt{a}\triangleq \min\left\{\frac{1}{2},\frac{2(1-\ln(2))}{M\sqrt{L_{\nabla f}}\diam(\mathcal{X})}\right\}$ and $\mathtt{b}\triangleq\frac{1-\ln(2)}{L_{\nabla f}\diam(\mathcal{X})^{2}}$. With the help of this lower bound, we are now able to establish the $\mathcal{O}(k^{-1})$ convergence rate in terms of the approximation error $h_{k}\triangleq f(x^{k})-f^{\ast}$. 

By convexity, we have $\gap(x^{k})\geq h_{k}$. Therefore, the lower bound for $\Delta_{k}$ in \eqref{eq:lowDelta} can be estimated as $\Delta_{k}\geq\min\{\mathtt{a} h_{k},\mathtt{b} h_{k}^{2}\}$. Hence,
\begin{equation}\label{eq:h}
h_{k+1}\leq h_{k}-\min\{\mathtt{a} h_{k},\mathtt{b} h_{k}^{2}\}\qquad\forall k\geq 0.
\end{equation}
Given this recursion, we can identify two phases characterizing the process $(h_{k})_{k\geq 0}$. In Phase I, the approximation error $h_k$ is at least $\mathtt{a}/\mathtt{b}$, and in Phase II the approximation error falls below this value. The cut-off value $\mathtt{a}/\mathtt{b}$ determines the nature of the recursion \eqref{eq:h}, and yields immediately an estimate for the iteration complexity of Variant 1 of Algorithm \ref{alg:AFW}.
\begin{theorem}\label{thm:1}
For all $\varepsilon>0$, define the stopping time 
\begin{equation}
N_{\varepsilon}(x^{0})\triangleq\inf\{k\geq 0\vert h_{k}\leq\varepsilon\}.
\end{equation}
Then, 
\begin{equation}\label{eq:N}
N_{\varepsilon}(x^{0})\leq \ceil[\bigg]{\frac{1}{\mathtt{a}}\ln\left(\frac{h_{0}\mathtt{b}}{\mathtt{a}}\right)} +\frac{L_{\nabla f}\diam(\mathcal{X})^{2}}{(1-\ln(2))\varepsilon}.
\end{equation}
\end{theorem}
The proof is in Section B.1 of the supplementary materials. 
\subsection{Analysis of Variant 2}
For the analysis of V2, we first need to establish well-posedness of the backtracking scheme $\mathtt{step}(f,v,x,g,\mathcal{L})$. Calling this routine at the position $x=x^{k}$ within the execution of Algorithm \ref{alg:AFW}, we require the search direction $v^{k}\triangleq s^{k}-x^{k}$, where $s^{k}=s(x^{k})$ is the target vector \eqref{eq:s}. Define 
\begin{equation}\label{eq:gamma}
\gamma_{k}\triangleq\max\{t\geq 0\vert x^{k}+t(s^{k}-x^{k})\in\mathcal{S}(x^{k})\}
\end{equation}
as the largest step size guaranteeing feasibility and decrease of the objective function value. Hence, for all $t\in[0,\gamma_{k}]$, we have $f(x^{k}+t(s^{k}-x^{k}))\leq f(x^{k})$, and $x^{k}+t(s^{k}-x^{k})\in\mathcal{X}$.

\begin{lemma}
Assume that $x^{k}\in\mathcal{S}(x^{0})$ for all $k\geq 0$. Then, for all $t\in[0,\gamma_{k}]$, it holds true that $x^{k}+t(s^{k}-x^{k})\in \mathcal{S}(x^{k})$, and 
\[
\norm{\nabla f(x^{k}+t(s^{k}-x^{k}))-\nabla f(x^{k})}_{2}\leq L_{\nabla f}t \norm{s^{k}-x^{k}}_{2}. 
\]
\end{lemma}
Using this estimate, a localized version of the celebrated descent Lemma reads as 
\begin{align*}
f(x^{k}+t(s^{k}-x^{k}))\leq Q(x^{k},t,L_{\nabla f}),
\end{align*}
for all $t\in[0,\gamma_{k}]$ for the quadratic model introduced in \eqref{eq:Q}. A direct minimization strategy for constructing a feasible step size policy based on the majorizing quadratic model $t\mapsto Q(x^{k},t,L_{\nabla f})$ would yield the step size 
\begin{equation}\label{eq:tau}
\tau_{k}(L_{\nabla f})\triangleq\min\left\{1,\frac{\gap(x^{k})}{L_{\nabla f}\norm{s^{k}-x^{k}}^{2}_{2}}\right\}.
\end{equation}
The main problem with this approach is that it needs the global parameter $L_{\nabla f}$. In practice this quantity is often hard to obtain and, furthermore, frequently numerically quite large. This in turn renders the step size to be inevitably small, leading to bad performance of the method. Our solution strategy is thus to implement the adaptive backtracking procedure calling the function $\mathtt{step}(f,v^{k},x^{k},\gap(x^{k}),\mathcal{L}_{k-1})$, which requires a local Lipschitz estimate $\mathcal{L}_{k-1}$, guaranteed to be smaller than $L_{\nabla f}$. This subroutine needs an initial estimate $\mathcal{L}_{-1}$ whose choice is explained in Section C of the supplementary materials.

Furthermore, Algorithm \ref{alg:backtrack} comes with two hyperparameters $\gamma_{d}<1<\gamma_{u}$. Practically efficient values for these parameters are reported to be $\gamma_{d}=0.9$ and $\gamma_{u}=2$ \cite{PedNegAskJag20}. Section C of the supplementary materials gives an upper bound estimate on the number of necessary function evaluations during a single execution of the backtracking procedure.
The $\mathcal{O}(k^{-1})$ convergence rate of V2 involves the \emph{dual objective function} 
$$
\psi(z)\triangleq -f^{\ast}(z)-H_{\mathcal{X}}(-z),
$$
where $f^{\ast}$ denotes the Fenchel conjugate of $f$ and $H_{\mathcal{X}}(c)\triangleq \sup_{x\in\mathcal{X}}\inner{x,c}$ is the support function over $\mathcal{X}$. Note that $\psi$ is concave and we have $f^{\ast}\triangleq \min_{x\in\mathcal{X}}f(x)=\max_{u}\psi(u)$. The following result is then similar to Theorem 3 in \citet{PedNegAskJag20}, and reports the $\mathcal{O}(k^{-1})$ convergence rate in terms of the sequence of approximation errors $h_{k}=f(x^{k})-f^{\ast}$.

\begin{theorem}\label{th:V2}
Let $f$ be an \ac{SC} function and $(x^{k})_{k\geq 0}$ generated by Variant 2 of Algorithm \ref{alg:AFW}.
Then
\begin{equation}
h_{k}\leq\frac{2\gap(x^{0})}{(k+1)(k+2)}+\frac{k\diam(\mathcal{X})^{2}}{(k+1)(k+2)}\bar{\mathcal{L}}_{k}
\end{equation}
where $\bar{\mathcal{L}}_{k}\triangleq\frac{1}{k}\sum_{i=0}^{k-1}\mathcal{L}_{i}$ is the arithmetic mean over all Lipschitz estimates computed by Algorithm \ref{alg:backtrack} during $k$ executions of the main protocol.
\end{theorem}
Since $\mathcal{L}_{k}\leq L_{\nabla f}$, this implies $h_{k}=\mathcal{O}(k^{-1})$. The proof of this Theorem can be found in supplementary material Section C.3.


\section{Inducing Linear Convergence}
\label{sec:fast}
In this section we use the construction of a \ac{LLOO} given in \citet{GarHaz16} to deduce an accelerated version of the base \ac{FW} Algorithm \ref{alg:AFW}. \citet{GarHaz16} give an explicit construction of such an oracle for \emph{any} polytope of form $\mathcal{X}=\{x\in\mathbb{R}^n\vert Ax= a,Bx\leq b\},$ where $A,B\in\mathbb{R}^{m\times n}$ and $a,b\in\mathbb{R}^{m}$. In the particularly important case where $\mathcal{X}$ is a standard unit simplex, the \ac{LLOO} can be constructed efficiently. In the supplementary material we describe the construction for this special geometry.
Let $\ball(x,r)\triangleq\{y\in\mathbb{R}^n\vert \norm{y-x}_{2}\leq r\}$ denote the Euclidean ball with radius $r>0$ and center $x$.

\begin{definition}[\citet{GarHaz16}, Def. 2.5]
A procedure $\mathcal{A}(x,r,c)$, where $x\in\mathcal{X},r>0,c\in\mathbb{R}^n,$ is a \ac{LLOO} with parameter $\rho\geq 1$ for the polytope $\mathcal{X}$ if $\mathcal{A}(x,r,c)$ returns a point $s\in\mathcal{X}$ such that for all $y\in\ball(x,r)\cap\mathcal{X}$
\begin{equation}
\inner{c,y}\geq\inner{c,s}\text{ and }\norm{x-s}_{2}\leq \rho r.
\end{equation}
\end{definition}

Assuming the availability of a procedure $\mathcal{A}(x,r,c)$ for any point $x\in\mathcal{X}$, we run Algorithm \ref{alg:LLOO}.
\begin{algorithm}
 \caption{LLOO-based convex optimization}
 \label{alg:LLOO}
 \begin{algorithmic}
\STATE {\bfseries Input: } $\mathcal{A}(x,r,c)$-LLOO with parameter $\rho\geq 1$ for polytope $\mathcal{X}$, $f\in\mathcal{F}_{M}$. $\sigma_{f}>0$ convexity parameter.\\
 $x^{0}\in\dom f\cap \mathcal{X}$, and let $h_{0}= f(x^{0})-f^{\ast}$, $c_{0}=1$.
 \FOR{k=0}
\STATE Obtain $r_{0}=\sqrt{\frac{6\gap(x^{0})}{\sigma_{f}}}$.
 \STATE Obtain $s^{0}=\mathcal{A}(x^{0},r_{0},\nabla f(x^{0}))$
 \STATE Set $\alpha_{0}=\frac{1}{1+\mathtt{e}^{0}}\min\{1,\frac{\gap(x^{0})}{\frac{4}{M^{2}}(\mathtt{e}^{0})^{2}}\}$, where $\mathtt{e}^{0}=\norm{s^{0}-x^{0}}_{x^{0}}$.
 \STATE Update $x^{1}=x^{0}+\alpha_{0}(s^{0}-x^{0})$ 
 \ENDFOR
\FOR{$k=1,\ldots$}
 \IF{$\gap(x^{k})>\varepsilon$}
\STATE Set $c_{k}=\exp\left(-\frac{1}{2}\sum_{i=0}^{k-1}\alpha_{i}\right)$
\STATE Set $\alpha_{k}=\min\{\frac{c_k\gap(x^0)}{\frac{4}{M^{2}}(\mathtt{e}^k)^2},1\}\frac{1}{1+\mathtt{e}^k}$
\STATE Set $r_{k}=r_{0}c_{k}$.
\STATE Obtain $s^{k}=\mathcal{A}(x^{k},r_{k},\nabla f(x^{k}))$
\STATE Set $x^{k+1}=x^{k}+\alpha_{k}(s^{k}-x^{k})$
\ENDIF
\ENDFOR
\end{algorithmic}
\end{algorithm}
Our analysis of this Algorithm departs from eq. \eqref{eq:up}, saying that for $\alpha\in(0,1/\mathtt{e}^{k})$ we have 
\begin{align*}
f(x^{k+1})\leq f(x^{k})+\alpha\inner{\nabla f(x^{k}),s^{k}-x^{k}}+\frac{4}{M^{2}}\omega_{\ast}(\alpha\mathtt{e}^{k}), 
\end{align*}
To exploit the power of the \ac{LLOO}, we need to control the right-hand side of the previous display by bounding carefully the radius of a single step made by the algorithm. Via a delicate induction argument, based on estimates for \ac{SC} functions and \citet{GarHaz16}, we obtain the announced linear convergence result.
\begin{theorem}
\label{thm:LCON}
Let $(x^{k})_{k\geq 0}$ be generated by Algorithm \ref{alg:LLOO}. Then for all $k\geq 0$ we have $x^{\ast}\in\ball(x^{k},r_{k})$ and 
\begin{equation}\label{eq:fast}
h_{k}\leq \gap(x^{0})\exp\left(-\frac{1}{2}\sum_{i=0}^{k-1}\alpha_{i}\right)
\end{equation}
In particular, defining $\bar{\alpha}\triangleq \min\{ \frac{\sigma_f}{6L_{\nabla f} \rho^2},1\}\frac{1}{1+ \sqrt{L_{\nabla f}}\frac{M\diam(\mathcal{X})}{2}}$, we see that $\alpha_{k}\geq\bar{\alpha}$, and therefore
\[
h_{k}\leq \gap(x^{0})\exp(-k\bar{\alpha}/2).
\]

\end{theorem}
The proof is in Section C of the supplementary materials. 




\section{Numerical Experiments}
\label{sec:numerics}
\begin{figure}[t]
\begin{center}
\centerline{\includegraphics[width=1.\columnwidth]{./Figures/icml.png}}
\caption{Performance of Algorithm \ref{alg:AFW} and \ref{alg:LLOO} in the Poisson inverse problem \eqref{eq:Poisson} (left column) and in the Portfolio optimization problem \eqref{eq:Portfolio} (right column). The top row reports the fraction of datasets for which a given relative error was achieved by an algorithm. The middle row presents the average ratio between number of iterations taken to reach a certain relative error to the method which reached this error in the minimal number of iterations. The bottom row presents the average time ratio between each method and that of the fastest method in order to reach a certain relative error.}
\label{fig:Portfolio}
\end{center}
\vskip -0.2in
\end{figure}

In the numerical experiments we tested the performance of Variant 1 (V1) and Variant 2 (V2) of Algorithm \ref{alg:AFW}, and compared them with the performance of Frank-Wolfe with {standard} step-size of $\frac{2}{k+2}$ (stand.), and step-size determined by exact {line-search} (Line-S.). As a further benchmark, the self-concordant Proximal-Newton (PN) of \citet{CevKyrTra15}, as implemented in the SCOPT package\footnote{\url{https://www.epfl.ch/labs/lions/technology/scopt/}}, is included. For the portfolio optimization problem, Algorithm \ref{alg:LLOO} is also implemented. All codes are written in Python 3, with packages for scientific computing NumPy 1.18.1 and SciPy 1.4.1. The experiments were conducted on a PC with Intel Core i5-7500 3.4GHzs, with a total of 16GB RAM. In both experiments the Frank-Wolfe based methods have been terminated after 50,000 iterations. Because of its higher computational complexity, we decided to stop PN after 1,000 iterations. Each algorithm was stopped early if the optimality gap in a given iteration was lower than $1e-10$. LLOO was only implemented for the portfolio selection problem, using the local linear oracle given in \cite{GarHaz16} for the simplex, as described in Section F in the supplementary materials.\footnote{The codes are accessible via \url{https://github.com/kamil-safin/SCFW }.}\\
For the portfolio optimization problem we used synthetic data, as in Section 6.4 of \cite{SunTran18}. The details of the data generating process are as follows. We generate matrix $R$ with given price ratios as: $R := ones(n, p) + N (0, 0.1)$, which allows the closing price to vary about 10\% between two consecutive periods. We used different sizes of matrix $R$: $(n, p) = (1000, 800), (1000, 1200)$, and $(1000, 1500)$ with 4 samples for each size. Hence, there are totally 12 datasets. The right column of Figure~\ref{fig:Portfolio} displays the performance of all Algorithms developed in this paper when averaged over all 12 samples. The detailed results for individual samples are reported in Section E of the supplementary materials. 

For the Poisson inverse problem we used the datasets a1a-a9a from the LIBSVM library \cite{LIBSVM}. Averaged performance of the \ac{FW} algorithms over these 9 datasets is displayed in the left column of Figure \ref{fig:Portfolio}. The detailed results for each individual data set are reported in Section E in the supplementary materials.

The comparison between all the algorithms is made by the construction of versions of performance profiles, following \citet{DolMor02}. In order to present the result, we first estimate $f^*$ by the maximal lower bound on the function value achieved by any of the algorithms, and compute the relative error attained by each of the methods at iteration $k$. More precisely, given the set of methods $\mathcal{S}$ and test problems $\mathcal{P}$, denote by $F_{ij}$ the function value attained by method $i$ on problem $j$. If $(x^{k}_{ij})_{k}$ denotes the sequence produced  by method $i$ on problem $j$, we define the relative error as 
$r^{k}_{ij}=\frac{f(x^{k}_{ij})-\min\{F_{sj}\vert s\in\mathcal{S}\}}{\min\{F_{sj}\vert s\in\mathcal{S}\}}$.
Now, for all methods $i\in\mathcal{S}$ and any relative error $\epsilon$, The  top row of Figure~\ref{fig:Portfolio} shows the proportion of datasets that achieves a relative error of $\epsilon$, that is $\rho_i(\epsilon)= \frac{1}{|\mathcal{P}|}\abs{\{j\in\mathcal{P}\vert \exists k,\; r^{k}_{ij}\leq \varepsilon\}}$.
We are also interested in comparing iteration complexity and CPU time between the methods. For that purpose, we define 
 $N_{ij}(\varepsilon)=\min\{k\geq 0\vert r^{k}_{ij}\leq \varepsilon\}$ and $T_{ij}(\varepsilon)$ as the first iteration and the minimal CPU time in which method $i\in\mathcal{S}$ to achieves a relative error $\varepsilon$ on problem $j\in\mathcal{P}$. 
In the second row of Figure~\ref{fig:Portfolio}, we use the following average iteration ratio $\tilde{\rho}_{i}(\varepsilon)=\frac{1}{|\mathcal{P}|}\sum_{j\in\mathcal{P}}\frac{N_{ij}(\varepsilon)}{\min\{N_{sj}(\varepsilon)\vert s\in\mathcal{S}\}}$ for comparing the iteration complexity of all the methods. 
In the third row of Figure~\ref{fig:Portfolio}, we display the average time ratio $\hat{\rho}_{i}(\varepsilon)=\frac{1}{|\mathcal{P}|}\sum_{j\in\mathcal{P}}\frac{T_{ij}(\varepsilon)}{\min\{T_{sj}(\varepsilon)\vert s\in\mathcal{S}\}}$ for comparing the computational time of all the methods. In both cases, as the average ratio is closer to 1 the performance of the method is closer to the best performance.

As expected, the top row of Figure~\ref{fig:Portfolio} show that PN obtains a lower relative error than any of the Frank-Wolfe based methods. However, for the Frank-Wolfe methods, none of the step-sizes choices has a clear advantage over the other in obtaining lower relative error values in both examples. Moreover, in the portfolio, while PN has the lowest iteration complexity, Variant 2 in fact achieves the best times for all values of relative error, followed by line-search, with LLOO matching the performance of the standard step-size. We remark that the standard step-size policy has no theoretical convergence guarantees, as the problem has no finite curvature constant (nor has a Lipschitz continuous gradient over the unit simplex).\\
In the Poisson inverse problem displayed in the left column of Figure \ref{fig:Portfolio}, Variant 2 obtains the best times for higher values of relative errors. The standard $2/(k+2)$ step size rule seems to be a good alternative here. We also implemented the convergent version with the step size $2/(k+3)$ derived in \citet{Odor16} and did not observe any significant difference in practice with the standard step-size policy. It is remarkable that the worst-case step-size policy V1 is competitive with the problem-specific $2/(k+3)$ policy. 


\section{Conclusion}
\ac{FW} is a much appraised \ac{FOM} for minimizing smooth convex functions over convex compact sets. The main merit of this method is the relative simplicity of its implementation and projection-free iterations. This yields great scalability properties making it a very attractive method for large-scale optimization. The price to pay for iteration-simplicity are, in general, slow convergence rates, and some sort of bounded curvature assumption.
In this work, we show that for \ac{SC} functions, which are neither strongly convex, nor have bounded curvature, we can obtain a novel step-size policy which, when coupled with local linear minimization oracles, features linear convergence rates. Under more standard assumptions of the feasible set given by \ac{LMO}, we provide two novel step-size policies which lead to standard sublinear convergence rate for minimizing general \ac{SC} functions by \ac{FW} algorithm.
In the future we plan to extend the results to the class of generalized self-concordant functions, as recently defined in \cite{SunTran18}. Other directions of interest for future research are inertial methods, and stochastic optimization. We will address these problems in the near future. 
An interesting avenue could be to incorporate away steps. This seems to be challenging since vanilla away step is not trivially to implement as the function might not be bounded and the away step may lead to a point outside of the domain of the objective.

\section*{Acknowledgements}
The authors sincerely thank Shoham Sabach for his contribution in the early stages of this project, including his part in developing the basic ideas used in this paper. We would also like to thank Quoc Tran-Dinh for sharing MATLAB codes. This research is supported by the COST Action CA16228 "European Network for Game Theory".



