# From Robot Self-Localization to Global-Localization:
An RSSI Based Approach.

Athanasios Lentzas Aristotle University
Thessaloniki, GreeceSchool of Informatics
Aristotle University of Thessaloniki
Thessaloniki, Greece<EMAIL_ADDRESS>Aristotle University
Thessaloniki, GreeceSchool of Informatics
Aristotle University of Thessaloniki
Thessaloniki, Greece Dimitris Vrakas Aristotle University
Thessaloniki, GreeceSchool of Informatics
Aristotle University of Thessaloniki
Thessaloniki, Greece<EMAIL_ADDRESS>
###### Abstract

Localization is a crucial task for autonomous mobile robots in order to
successfully move to goal locations in their environment. Usually, this is
done in a robot-centric manner, where the robot maintains a map with its body
in the center. In swarm robotics applications, where a group of robots needs
to coordinate in order to achieve their common goals, robot-centric
localization will not suffice as each member of the swarm has its own frame of
reference. One way to deal with this problem is to create, maintain and share
a common map (global coordinate system), among the members of the swarm. This
paper presents an approach to global localization for a group of robots in
unknown, GPS and landmark free environments. The main idea relies on members
of the swarm staying still and acting as beacons, emitting electromagnetic
signals. These stationary robots form a global frame of reference and the rest
of the group localize themselves in it using the Received Signal Strength
Indicator (RSSI). The proposed method is evaluated, and the results obtained
from the experiments are promising.

##  1 Introduction

Robotic Swarms is a topic with significant interest, especially with the
increased use of large groups of unmanned air or ground vehicles [14]. Self
localization and navigation is a crucial task for autonomous robots [4], both
on indoor and outdoor navigation. Depending on the application of the robotic
swarm, the robots may need the ability to discover their location.

Non robotic applications can also take advantage of indoor localization
techniques. In case of emergency evacuation, especially on crowded buildings
like malls, airports etc., it is beneficial to know the location of people
inside the building as this can help with faster and safer evacuation. In
order to achieve that, specific devices or smartphones can be used that take
advantage of the proposed localization scheme.

Localization is the ability of a robot to know its position and orientation.
The location could be relative to other robots or absolute on a common
coordinate system. While operating outdoors, a robot can rely on GPS to
localize itself. The aforementioned sensor can not be used on indoor
environments. Simultaneous localization and mapping (SLAM) algorithms are
extensively used [3], although when small and relatively simple robots are
needed, their high complexity and computational cost is a major drawback.
Several techniques exploiting beacons and landmarks for localization have been
proposed in the literature as well as techniques using inertial navigation
systems, magnetic, sound and optic based navigation [15].

Inertial localization is also a common approach. Based on odometers,
accelerometers and gyroscopes, one can determine the orientation and direction
of the robot [8]. Although the results provided are robust against
environmental changes, this technique is prone to error accumulation [16].
Kalman filters have been expensively used in order to improve accuracy [9].

Radio frequency localization is the approach usually preferred. Beacons can
cover a wide area, radio waves can penetrate most materials while the
installation cost is relatively low [12]. Received Signal Strength Indicator
(RSSI) is mainly exploited for robot localization [20, 28, 2]. Time of arrival
and time difference of arrival of two signals that are known to have different
propagation is also a common approach for localization.

Localization based on RSSI is proposed in the Ladybug algorithm [11]. A robot
equipped with sensors able to measure the strength of the received
electromagnetic signal, is able to identify the location of the source of the
signal and navigate to it. The source could be either a beacon or another
robot. The LadyBug Algorithm was effective and had numerous benefits compared
to similar approaches, such as I-Bug [23]. RSSI could be an efficient solution
when deploying a robotic swarm on GPS denying environment. Our motivation was
to propose an approach where a robotic swarm could be able to extract the
location of each robot in the same coordinate system using local sensing only,
allowing the LadyBug algorithm to be implemented on the swarm. Sharing a
common coordinate system is crucial for tasks such as self-assembly.

The benefits of localizing the swarm in a common coordinate system are
presented in [19]. The ability of the swarm to self-assemble on a specific
shape was realized by placing four robots on a specific layout, marking the
coordinate system and using trillateration. While the proposed solution is
efficient, it’s main drawback is the limitation to two dimensions and the
requirement of hand placing the four initial robots. Our incentive was to
propose an approach where the robots will use local sensing but will have the
ability to localize in a global system without the aforementioned limitations.

Our approach would have numerous benefits such as: a) the swarm will share a
common localization scheme, sharing the same map. b) The localization scheme
is based on three members of the swarm equipped with a beacon instead of four.
c) No manual placement is required. The beacons could be anywhere on the
operation area. d) Our approach is able to localize in three dimensions. e)
The proposed solution is robust as, in case of a possible failure of a beacon,
another robot from the swarm can replace it.

##  2 Related Work

Radio frequency identification technology (RFID) is used in [24] to localize
robots navigating indoors. The passive RFID tags installed, divided the area
into a grid. As a robot, equipped with an RFID reader, explores the area, the
reader reads the tag. The position is then estimated by correlating the ID of
the tag with a map containing the localization information of the tags.
Despite the high accuracy, the need for a map containing the ID and location
of the tags limits the practical applications.

Pseudolites (i.e. pseudo satellites) have been proposed for indoor
localization [25]. The main idea is to receive the GPS signal and transmit it
indoors using signal repeaters [27]. The major drawback of this approach is
the increased cost of the network installation. Additionally, in case the main
signal receiver fails the whole network is unusable.

A localization method based on RSSI of heterogeneous sources (i.e. WLAN, GSM
etc) is presented in [22]. By analyzing the fingerprint and strength of the
received signal, the robot is able to localize itself by comparing it with a
fingerprint map. While exploiting already existing infrastructures, negating
the need for a network deployment, the main disadvantage is the need for a
fingerprint map covering all the area of interest. RFID sensors can also be
used for indoor localization. In [7], authors presented a simulator that
allows modeling environments and testing the deployment of RFID solutions.
Tags and antennas are placed and RSSI is exploited to find the location of
each tag.

A distributed localization method based on swarm intelligence algorithms is
presented in [5]. Particle swarm optimization and an approach based on
backtracking search algorithm is proposed. The proposed method operates in
three stages. The first stage allows a robot to estimate the distance from a
reference node using the Sum-Dist algorithm [10]. The second stage estimates
the position of the robot using the min-max technique[21, 10]. Finally the
third stage the accuracy improves by re-evaluating the position based on the
position of the neighbors.

##  3 Problem Statement

The localization of the robot is based on RSSI. Given an electromagnetic
signal, the coordinates of the source in 3D space can be calculated in the
robocentric system and the in the global coordinate system. In order to
achieve our goal, three beacons are employed, creating a 3D common coordinate
system for all robots. In our work, the following assumptions were made: a)
All robots are equipped with sensors capable of reading the RSSI of a received
signal. b) Two robots act as anchors. c) A beacon/robot is placed at the
center of the coordinate system.

Figure 1: The position of the four sensors on the robot.

The main localization scheme employs three receivers and one transceiver. The
transceiver allows the robot to act as a beacon or exchange information with
the rest of the swarm. The receivers are placed on the front, right and left
of a robot with radius r and the transceiver on the center, as seen in Fig. 1.
Each sensor is able to identify both the source of the signal and the signal
strength. The signal strength (S) is inversely proportional to the square of
the distance (d) between the sensor and the source (1). Using (1) the distance
between the source and each sensor can be calculated (1).

$\displaystyle S$ $\displaystyle=\frac{1}{d^{2}}$ $\displaystyle d$
$\displaystyle=\sqrt{\frac{1}{S}}$ (1)

The coordinates of each sensor, in the robot-centric system, can be seen on
Table 1. Let $(x_{s},y_{s},z_{s})$ be the coordinates of the signal source in
the robot-centric system. Using (1) the distance between each sensor (center,
front, left, right) is calculated.

Sensor | Robocentric Coordinates (x,y) | Enchanced Localization Coordinates
---|---|---
$S_{0}$ | (0,0,0) | (0,0)
$S_{L}$ | (0,r,0) | (0,-r)
$S_{R}$ | (0,-r,0) | (0,r)
$S_{F}$ | (r,0,0) | (r,0)
$S_{B}$ | - | (-r,0)
Table 1: Position of the sensors

$\displaystyle d_{0}^{2}$ $\displaystyle=x_{s}^{2}+y_{s}^{2}+z_{s}^{2}$
$\displaystyle d_{f}^{2}$
$\displaystyle=x_{s}^{2}+2rx+r^{2}+y_{s}^{2}+z_{s}^{2}$ $\displaystyle
d_{l}^{2}$ $\displaystyle=x_{s}^{2}+y_{s}^{2}+z_{s}^{2}+2ry+r^{2}$
$\displaystyle d_{r}^{2}$
$\displaystyle=x_{s}^{2}+y_{s}^{2}+z_{s}^{2}-2ry+r^{2}$ (2)

Solving the equation system (3), the coordinates of the source, $x_{s}$ and
$y_{s}$ respectively, can be calculated. By replacing the values of $x_{s}$
and $y_{s}$ on (3) the $z_{s}$ coordinate can be calculated. Similarly (and
knowing the coordinates of the source) the coordinates of the beacon
$(x_{b},y_{b},z_{b})$ and the anchor robots $(x_{a1},y_{a1},z_{a1})$,
$(x_{a2},y_{a2},z_{a2})$ can be calculated.

###  3.1 Enhanced Localization

In order to further enhance the localization procedure, redundant sensors were
added. Our implementation uses five sensors. All the receivers are placed on
the perimeter of the robot thus creating four sensor groups (triangle formed
by the transceiver in the center and two receivers on the perimeter). The
updated robot-centric coordinates can be seen on Table 1. The estimates of the
four groups are averaged, resulting in more accurate estimation of the
source’s position.

###  3.2 Global Localization

As already described in the previous section, the robot is able to calculate
the coordinates of the three beacons in the robocentric system. Knowing their
relative positions, the global position of the robot can be calculated. Let A
be the center of the global coordinate system and B, C the two beacons that
define axis X and Y respectively.

Using vectors $\vec{AB}$ and $\vec{AC}$ we define three new vectors using the
outer products, as seen in (3), and normalize them.

$\displaystyle\vec{z}$ $\displaystyle=\vec{AB}\otimes\vec{AC}$
$\displaystyle\vec{y}$ $\displaystyle=\vec{z}\otimes\vec{AB}$
$\displaystyle\vec{x}$ $\displaystyle=\vec{y}\otimes\vec{z}$ (3)

Let $\vec{r}$ be the vector between the center of the robot R and A. We can
now calculate the coordinates of the robot in the global system defined by
$\vec{x}$, $\vec{y}$ and $\vec{z}$. The coordinates of the robot,
$(x_{r},y_{r},z_{r})$, are the dot product of $\vec{r}$ and $\vec{x}$,
$\vec{y}$, $\vec{z}$ respectively.

##  4 Experimental Results

In order to evaluate the performance of our localization approach we
implemented the algorithm in the Webots Open Source Robot simulator [26]. A
swarm of four identical robots was used, where three of them were serving as
beacons. A generic round shaped robot with 20cm radius was employed, equipped
with four radio receivers capable of identifying the signal strength on its
perimeter and a radio transceiver on the center.

The three robots acting as beacons were constantly emitting a radio signal
with a unique ID, allowing the robot to identify the source of the signal. In
order to further increase the accuracy, each beacon was transmitting on a
different channel reducing the generated noise. For each source the RSSI of
100 packets was averaged reducing the noise generated error. The robot was
cycling through the three predefined channels allowing it to identify and
locate the robocentric position of each beacon. Lastly the global coordinates
were calculated and reported.

During the experiment, the robot was static while calculating it’s position.
After a successfully calculation it was moving randomly and the new position
was calculated and reported. This process was repeated 10 times. The error
(i.e. euclidean distance between the real and calculated global location) for
each location was averaged and used as a metric. In a noiseless environment,
the calculated position was the same with the real one, no matter the distance
from the beacons or the orientation of the robot.

In order to evaluate our approach in more real-like conditions, noise was
gradually introduced into the experiment. In total the process described above
was repeated 5 times with 5 different noise levels: 10%, 20%, 30%, 40% and
50%. Noise strength is the standard deviation of the Gaussian noise added to
the signal strength. The rest of the conditions were controlled (i.e. the
random seed) making sure the experiment would be the same, providing
comparable results. As seen in Fig. 2a, for 10% noise the error was 0.6m while
for 20% noise the calculated position was 1.4m off. Further increasing the
noise had a bigger impact.

(a) Error per noise level

(b) Error per distance from the beacons

Figure 2: Localization error.

As noise is related to the distance of the source, an experiment was performed
where the distance between the robot and beacons was gradually increased. The
three beacons were hand placed forming an isosceles triangle. This placement
guaranteed that the distance between the robot and two beacons would be the
same. The side of the triangle formed was 4m and the noise level was set to
10%. The experiment was repeated 10 times and the results were averaged. As
seen in Fig. 2b the error was 0.1m when the robot was closer to the source.
The best results were observed when the robot was between 3m and 3.5m from the
beacons. Increasing the distance between the robot and source results in
increased noise level, making it harder to correctly identify the position of
the robot. It is worth mentioning that for distances up to almost 6m the
calculated position error was less than 1m.

An experiment was also performed using a small swarm of three robots. The
first robot was placed away from the beacons, the second close to the two of
them while the third close to the three beacon robots. The robots reported
their global position 10 times with a time delay added between each
calculation. The noise level was set to 10%. Further noise was introduced, as
each robot was emitting a random signal, simulating an environment where the
robots communicate with each other. In Fig. 3 the position of each robot can
be seen. The results show that even with the extra noise (i.e. robots
communicating with each other), the localization process is not heavily
affected. It is worth mentioning that the distance between the robot and the
beacons has an impact on the calculations (as already discussed).

Figure 3: Calculated position of each robot.

##  5 Conclusion & Future Work

In this paper, a novel global localization method was presented. Our approach
exploits RSSI from a signal emitted by three beacons, to localize a robot on
the coordinate system formed by those beacons.

Our approach was evaluated in a noiseless environment as well as in an
environment with different noise levels. While for lower noise levels the
results were accepted, for increased noise levels the error was bigger. This
requires further experiments as noise levels in the real world are affected by
many factors, such as the materials of the obstacles, reflectance etc.

Additionally, the impact of the distance between the robot and the source of
the signal was investigated. Placing the robot close to the beacons negatively
affects the ability to localize itself. Increasing the distance up to a
certain threshold has low impact. Further increasing the distance produced
increased error.

Although a small swarm was simulated, allowing the observation of how multiple
robots communicating will impact the performance, the robots were not
exchanging any location related information. In the near future our focus will
be on investigating ways that will allow the robots to correlate their
individual localization information and create a more accurate global map.

Furthermore, the current experiments were performed in a simulated
environment. Even with noise added it is not as close as the real world.
Experiments will be performed on a real robot. That evaluation is crucial as
real world noise can not be reproduced easily (i.e. reflection of the signal,
non-uniform noise etc.).

Another important part of our future work is to integrate the localization
presented with the LadyBug algorithm [23]. As LadyBug provides robust results,
the proposed method will allow its application to robotic swarms.

##  6 Acknowledgements

This research was carried out as part of the project «Beacon: An Intelligent
Dynamic Signage System for Emergency Situations» (Project code: KMP6-0094837)
under the framework of the Action «Investment Plans of Innovation» of the
Operational Program «Central Macedonia 2014 2020», which is co-funded by the
European Regional Development Fund and Greece.

## References

  * [1]
  * [2] Paolo Barsocchi, Stefano Lenzi, Stefano Chessa & Gaetano Giunta (2009): _A Novel Approach to Indoor RSSI Localization by Automatic Calibration of the Wireless Propagation Model_.  In: VTC Spring 2009 - IEEE 69th Vehicular Technology Conference, pp. 1–5, 10.1109/VETECS.2009.5073315.
  * [3] G. Bresson, Z. Alsayed, L. Yu & S. Glaser (2017): _Simultaneous Localization and Mapping: A Survey of Current Trends in Autonomous Driving_.  IEEE Transactions on Intelligent Vehicles 2(3), pp. 194–220, 10.1109/TIV.2017.2749181.
  * [4] I.J. Cox (1991): _Blanche-an experiment in guidance and navigation of an autonomous robot vehicle_.  IEEE Transactions on Robotics and Automation 7(2), pp. 193–204, 10.1109/70.75902.
  * [5] Alan Oliveira de Sa, Nadia Nedjah & Luiza de Macedo Mourelle (2016): _Distributed efficient localization in swarm robotic systems using swarm intelligence algorithms_.  Neurocomputing 172, pp. 322–336, 10.1016/j.neucom.2015.03.099.  Available at https://www.sciencedirect.com/science/article/pii/S0925231215010498.
  * [6] Qian Dong & Waltenegus Dargie (2012): _Evaluation of the reliability of RSSI for indoor localization_.  In: 2012 International Conference on Wireless Communications in Underground and Confined Areas, pp. 1–6, 10.1109/ICWCUCA.2012.6402492.
  * [7] Salvatore D’Avella, Matteo Unetti & Paolo Tripicchio (2022): _RFID Gazebo-Based Simulator With RSSI and Phase Signals for UHF Tags Localization and Tracking_.  IEEE Access 10, pp. 22150–22160, 10.1109/ACCESS.2022.3152199.
  * [8] Giovanni Fusco & James M. Coughlan (2018): _Indoor Localization Using Computer Vision and Visual-Inertial Odometry_.  In Klaus Miesenberger & Georgios Kouroupetroglou, editors: Computers Helping People with Special Needs, Springer International Publishing, Cham, pp. 86–93, 10.1007/978-3-319-94274-2_13.
  * [9] Guanghui Hu, Weizhi Zhang, Hong Wan & Xinxin Li (2020): _Improving the Heading Accuracy in Indoor Pedestrian Navigation Based on a Decision Tree and Kalman Filter_.  Sensors 20(6), 10.3390/s20061578.  Available at https://www.mdpi.com/1424-8220/20/6/1578.
  * [10] Koen Langendoen & Niels Reijers (2003): _Distributed localization in wireless sensor networks: a quantitative comparison_.  Computer Networks 43(4), pp. 499–518, 10.1016/S1389-1286(03)00356-6.  Available at https://www.sciencedirect.com/science/article/pii/S1389128603003566.  Wireless Sensor Networks.
  * [11] Athanasios Lentzas & Dimitris Vrakas (2020): _LadyBug. An Intensity based Localization Bug Algorithm_.  In: 2020 25th IEEE International Conference on Emerging Technologies and Factory Automation (ETFA), 1, pp. 682–689, 10.1109/ETFA46521.2020.9212115.
  * [12] Junjie Liu (2014): _Survey of Wireless Based Indoor LocalizationTechnologies_.  Department of Science and Engineering Washington University.
  * [13] O. Michel (2004): _Webots: Professional Mobile Robot Simulation_.  Journal of Advanced Robotics Systems 1(1), pp. 39–42.  Available at http://www.ars-journal.com/International-Journal-of-Advanced-Robotic-Systems/Volume-1/39-42.pdf.
  * [14] Nadia Nedjah & Luneque Silva Junior (2019): _Review of methodologies and tasks in swarm robotics towards standardization_.  Swarm and Evolutionary Computation 50, p. 100565, 10.1016/j.swevo.2019.100565.  Available at https://www.sciencedirect.com/science/article/pii/S2210650217308398.
  * [15] Huthaifa Obeidat, Wafa Shuaieb, Omar Obeidat & Raed Abd-Alhameed (2021): _A Review of Indoor Localization Techniques and Wireless Technologies_.  Wireless Personal Communications, 10.1007/s11277-021-08209-5.
  * [16] Enrico Petritoli, Fabio Leccese & Mariagrazia Leccisi (2019): _Inertial Navigation Systems for UAV: Uncertainty and Error Measurements_.  In: 2019 IEEE 5th International Workshop on Metrology for AeroSpace (MetroAeroSpace), pp. 1–5, 10.1109/MetroAeroSpace.2019.8869618.
  * [17] Joseph A. Rothermich, M. İhsan Ecemiş & Paolo Gaudiano (2005): _Distributed Localization and Mapping with a Robotic Swarm_.  In Erol Şahin & William M. Spears, editors: Swarm Robotics, Springer Berlin Heidelberg, Berlin, Heidelberg, pp. 58–69, 10.1007/978-3-540-30552-1_6.
  * [18] Stergios I. Roumeliotis & George A. Bekey (2000): _Distributed Multi-Robot Localization_ , pp. 179–188.  Springer Japan, Tokyo, 10.1007/978-4-431-67919-6_17.
  * [19] Michael Rubenstein, Alejandro Cornejo & Radhika Nagpal (2014): _Programmable self-assembly in a thousand-robot swarm_.  Science 345(6198), pp. 795–799, 10.1126/science.1254295.
  * [20] Jirapat Sangthong, Jutamas Thongkam & Sathapom Promwong (2020): _Indoor Wireless Sensor Network Localization Using RSSI Based Weighting Algorithm Method_.  In: 2020 6th International Conference on Engineering, Applied Sciences and Technology (ICEAST), pp. 1–4, 10.1109/ICEAST50382.2020.9165300.
  * [21] Andreas Savvides, Heemin Park & Mani B. Srivastava (2002): _The Bits and Flops of the N-Hop Multilateration Primitive for Node Localization Problems_.  In: Proceedings of the 1st ACM International Workshop on Wireless Sensor Networks and Applications, WSNA ’02, Association for Computing Machinery, New York, NY, USA, p. 112–121, 10.1145/570738.570755.
  * [22] M. Stella, M. Russo & D. Begušić (2014): _Fingerprinting based localization in heterogeneous wireless networks_.  Expert Systems with Applications 41(15), pp. 6738–6747, 10.1016/j.eswa.2014.05.016.  Available at https://www.sciencedirect.com/science/article/pii/S0957417414002966.
  * [23] Kamilah Taylor & Steven M. LaValle (2009): _I-Bug: An intensity-based bug algorithm_.  In: 2009 IEEE International Conference on Robotics and Automation, pp. 3981–3986, 10.1109/ROBOT.2009.5152728.
  * [24] R. Tesoriero, R. Tebar, J.A. Gallud, M.D. Lozano & V.M.R. Penichet (2010): _Improving location awareness in indoor spaces using RFID technology_.  Expert Systems with Applications 37(1), pp. 894–898, 10.1016/j.eswa.2009.05.062.  Available at https://www.sciencedirect.com/science/article/pii/S095741740900503X.
  * [25] Xiaoguang Wan & Xingqun Zhan (2011): _The Research of Indoor Navigation System using Pseudolites_.  Procedia Engineering 15, pp. 1446–1450, 10.1016/j.proeng.2011.08.268.  Available at https://www.sciencedirect.com/science/article/pii/S1877705811017693.  CEIS 2011.
  * [26] Webots: _http://www.cyberbotics.com_.  Available at http://www.cyberbotics.com.  Commercial Mobile Robot Simulation Software.
  * [27] Rui Xu, Wu Chen, Ying Xu & Shengyue Ji (2015): _A New Indoor Positioning System Architecture Using GPS Signals_.  Sensors 15(5), pp. 10074–10087, 10.3390/s150510074.  Available at https://www.mdpi.com/1424-8220/15/5/10074.
  * [28] Zheng Yang, Zimu Zhou & Yunhao Liu (2013): _From RSSI to CSI: Indoor Localization via Channel Response_.  ACM Comput. Surv. 46(2), 10.1145/2543581.2543592.


# Grad-GradaGrad? A Non-Monotone Adaptive Stochastic Gradient Method

Aaron Defazio
Meta AI
Baoyu Zhou
Lehigh University
Lin Xiao
Meta AI

###### Abstract

The classical AdaGrad method adapts the learning rate by dividing by the
square root of a sum of squared gradients. Because this sum on the denominator
is increasing, the method can only decrease step sizes over time, and requires
a learning rate scaling hyper-parameter to be carefully tuned. To overcome
this restriction, we introduce GradaGrad, a method in the same family that
naturally grows or shrinks the learning rate based on a different accumulation
in the denominator, one that can both increase and decrease. We show that it
obeys a similar convergence rate as AdaGrad and demonstrate its non-monotone
adaptation capability with experiments.

##  1 Introduction

We consider stochastic optimization problems of the form

$\operatorname*{minimize}_{x\in\mathbb{R}^{d}}\leavevmode\nobreak\
f(x):=\mathbb{E}_{\xi}[f(x,\xi)],$

where $\xi$ is a random variable. We assume that each $f(x,\cdot)$ is convex
in $x$ but potentially non-smooth. Problems of this class are most often
addressed with stochastic gradient methods (SGD). Aiming to accelerate the
convergence of stochastic gradient methods in practice, adaptive rules for
adjusting the learning rate have a long history (e.g., Kesten, 1958) and have
been extensively studied more recently. The stochastic gradient descent
method, a member of the Mirror Descent (MD) family takes the general form:

$x_{k+1}=x_{k}-A_{k+1}^{-1}g_{k},$ (1)

where $A_{k+1}$ is a dynamically adjusted positive definite matrix. Due to the
high dimensionality of modern machine learning applications, $A_{k+1}$ is
usually taken to be a diagonal matrix. AdaGrad (Duchi et al., 2011) uses the
choice:

$A_{i,k+1}^{-1}=\frac{\gamma}{\sqrt{\textstyle\sum_{t=0}^{k}g_{i,t}^{2}}},$
(2)

where $A_{i,k+1}$ is the $i$th diagonal entry of $A_{k+1}$, $g_{i,t}$ is the
$i$th coordinate of the stochastic gradient $g_{t}$, and $\gamma$ is a hyper-
parameter. Another popular and effective algorithm is Adam (Kingma & Ba,
2014), which uses an exponential weighted average instead of a summation, in
combination with momentum. Both methods have seen widespread adoption and have
been the focus of extensive followup work (e.g., Tieleman & Hinton, 2012;
Reddi et al., 2019; Li & Orabona, 2019).

A fundamental limitation of AdaGrad is _monotone adaptation_. In other words,
since the sum in (2) is increasing, it can only decrease the learning rate
over time. Therefore, the global learning rate parameter $\gamma$ must be
carefully tuned. In addition, once the learning rate becomes very small, it
cannot pick up speed again even if the landscape becomes flat later on.

Input: $\gamma_{0}$, $\rho$, $r$, $x_{0}$

Initialize $g_{-1}=0$, $\alpha_{0}=0$

for $k=0$ to $n$ do

$g_{k}=\nabla f(x_{k},\xi_{k})$

$v_{k}=\left\|g_{k}\right\|^{2}-\rho\left\langle g_{k},g_{k-1}\right\rangle$

if $v_{k}\geq 0$ then

$\gamma_{k+1}=\gamma_{k}$

$\alpha_{k+1}=\alpha_{k}+v_{k}$

else

$v_{k}=\max\left(v_{k},-r\alpha_{k}\right)$

$\gamma_{k+1}=\gamma_{k}\sqrt{1-v_{k}/\alpha_{k}}$

$\alpha_{k+1}=\alpha_{k}$

end if

$A_{k+1}=\sqrt{\alpha_{k+1}}/\gamma_{k+1}$

$x_{k+1}=x_{k}-A_{k+1}^{-1}g_{k}$

end for

algocfGradaGrad

(Simplified Scalar variant)

figureThe GradaGrad update (Equation 9) results in a smooth scaling of the
learning rate when negative values of $v_{k}$ are encountered, compared to the
naive solution from Equation (3) which quickly explodes.

In this paper, we propose a _non-monotone_ adaptive stochastic gradient method
called Grad-GradaGrad, or GradaGrad for short. The basic idea is to append an
inner product term $\langle g_{t},g_{t-1}\rangle$ to the summands in (2).
Specifically,

$A_{i,k+1}\propto\sqrt{\textstyle\sum_{t=0}^{k}\left(g_{i,t}^{2}-\rho
g_{i,t}g_{i,t-1}\right)},$ (3)

where $\rho>0$ is a constant that can be tuned. The expectation of
$-g_{i,t}g_{i,t-1}$ is the hyper-gradient of $f$ with respect to the learning
rate for the $i$th coordinate (see Section 2). Intuitively, if consecutive
stochastic gradients are highly positively correlated, then $g_{i,t}^{2}-\rho
g_{i,t}g_{i,t-1}$ can be negative and thus the learning rate will be
increased. Similarly, if consecutive gradients are negatively correlated, then
the learning rate will be decreased. In the “goldilocks” zone, the hyper-
gradient term wobbles around 0 and the learning rate changes in a similar
manner to AdaGrad. Compared to AdaGrad, our method is able to rapidly increase
the learning rate when it is too small (Figure 1).

However, applying the idea above directly runs into several problems. The most
apparent problem is that the sum under the square root may become a negative
number. A more subtle problem is how to ensure that the sequences of
$A_{i,k+1}$ increase, in a non-monotone manner, at an appropriate rate to
guarantee convergence. After laying out the basic ideas in Section 3, we show
in Section 4 how the GradaGrad method addresses these problems with
reparametrization and other techniques. As a preview, Algorithm 1 shows a
simplified version of GradaGrad using a scalar learning rate.

In the rest of Section 4, we show that GradaGrad obeys a similar convergence
bound as the AdaGrad method. In Section 5, preliminary experiments demonstrate
that our method matches the practical performance of several benchmark methods
and has the advantage of non-monotone adaptation.

#### Notation and Assumptions

Throughout this paper, we let $x_{*}$ denote any minimizer of $f$ and $x_{0}$
be the initial point. We use the convention that $g_{-1}$ is the zero vector.
The notation $B\succ 0$ means that $B$ is a symmetric and positive definite
matrix. For $B\succ 0$, we define $\left\|x\right\|_{B}=\sqrt{x^{T}Bx}$ and
$\left\langle x,y\right\rangle_{B}=x^{T}By$. Some of our results will use the
Lipschitz smoothness of $f$, in which case we denote the smoothness constant
as $L$. We define the filtration
$\Phi=\left(\mathcal{F}_{1},\mathcal{F}_{2},\dots\right)$, where
$\mathcal{F}_{k}$ is the $\sigma$-algebra generated by a sequence of random
variables $\\{\xi_{0},\ldots,\xi_{k-1}\\}$. The notation
$\mathbb{E}_{\xi_{k}}[\,\cdot\,]$ and
$\mathbb{E}[\,\cdot\,|\,\mathcal{F}_{k}]$ denote the expectation conditioned
on $\mathcal{F}_{k}$.

(a) AdaGrad (b) GradaGrad

Figure 1:  (a) Monotone learning rates cause AdaGrad to converge slowly when
starting with a small initial learning rate; (b) GradaGrad automatically
adapts to a poor initial learning rate in a non-monotone manner, leading to
faster convergence.

##  2 Related Work

Adaptive stochastic gradient methods have seen heavy investigation on a number
of fronts, including non-monotone methods that can both increase and decrease
the learning rate. Schaul et al. (2013) proposed a method that combines the
estimation of gradient variance and local curvature. Orabona & Tommasi (2017)
proposed a method to adjust the learning rate based on a “coin-betting”
strategy. Methods that maintain multiple learning rates, such as Metagrad (van
Erven & Koolen, 2016) are also potentially non-monotonic. Another line of work
that can both increase and decrease the learning rate are stochastic line
search methods (Vaswani et al., 2019; Paquette & Scheinberg, 2020; Zhang et
al., 2020). None of the existing strategies can be considered “drop-in”
replacements for AdaGrad in the same way our approach is, as they require
either additional gradient evaluations, additional memory overhead or
knowledge of problem-dependent constants that make them harder to use in
practice than the AdaGrad family of methods.

Methods in the Polyak class (Loizou et al., 2021) compute learning rates using
the function value sub-optimality $f(x_{k})-f(x_{*})$, or in the stochastic
case $f(x_{k},\xi_{k})-f^{*}_{\xi_{k}}$ where $f^{*}_{\xi_{k}}$ is the
minimial value of $f(\cdot,\xi_{k})$. These methods can exhibit impressive
performance but rely on knowledge or estimates of the minimal function value
$f_{*}$ or $f^{*}_{\xi_{k}}$, which is not available for most problems.

The idea of adjusting the learning rate based on inner products of consecutive
stochastic gradients goes back to Kesten (1958), and it was extended to the
multi-dimensional case by Delyon & Juditsky (1993). Mirzoakhmedov & Uryasev
(1983) showed that such adaptation schemes can be interpreted as stochastic
hyper-gradient methods. Specifically, suppose $d_{k}$ is a search direction
depending on the random variable $\xi_{k}$ (e.g., $g_{k}$ itself or combined
with momentum). Let’s define a merit function of the learning rate at each
iteration

$\phi_{k}(\eta)=\mathbb{E}\left[f(x_{k+1})\,|\,\mathcal{F}_{k}\right]=\mathbb{E}\left[f(x_{k}-\eta
d_{k})\,|\,\mathcal{F}_{k}\right].$

The hypergradient of $\phi_{k}$ with respect to $\eta$ can be derived as

$\displaystyle\frac{\partial\phi_{k}}{\partial\eta}$
$\displaystyle=\mathbb{E}\left[\frac{\partial}{\partial\eta}f(x_{k}-\eta
d_{k})\,|\,\mathcal{F}_{k}\right]=\mathbb{E}\left[\langle\nabla f(x_{k}-\eta
d_{k}),-d_{k}\rangle\,|\,\mathcal{F}_{k}\right]$
$\displaystyle=\mathbb{E}\left[-\langle
g_{k+1},d_{k}\rangle\,|\,\mathcal{F}_{k}\right],$ (4)

where we assume necessary technical conditions for changing the order of
differentiation and expectation. Therefore, $-\langle g_{k+1},d_{k}\rangle$
can be viewed as a stochastic hypergradient. If $d_{k}=g_{k}$, then it becomes
$-\langle g_{k+1},g_{k}\rangle$. If in addition the learning rate is different
for each coordinate, say $\eta_{i}$ for the $i$th coordinate, then it can be
shown that

$\frac{\partial\phi_{k}}{\partial\eta_{i}}=\mathbb{E}\left[-\langle
g_{i,k+1},g_{i,k}\rangle\,|\,\mathcal{F}_{k}\right].$

These are exactly the terms, after multiplying by a constant $\rho$, that are
added to the summands under the square root in (3), although the time-step is
shifted backwards by one. In the Adam-style variant of GradaGrad we develop,
$d_{k}$ is the convex combination of $g_{k}$ and a momentum term, and the
additional terms in calculating $A_{i,k}$ becomes $-\rho g_{i,t}d_{i,t-1}$.

Adaptive stochastic gradient methods based on hyper-gradient descent have been
studied in the machine learning literature before, including Jacobs (1988),
Sutton (1992), Schraudolph (1999), Mahmood et al. (2012), and more recently by
Baydin et al. (2018). The major difference between our work and these previous
work is that, although our adaptation scheme carries the intuition of hyper-
gradient, our method is not really a hyper-gradient method and our analysis
does not rely on the hyper-gradient interpretation.

The “without descent” adaptive method of Malitsky & Mishchenko (2019) is the
most direct inspiration for our work. In a footnote, they suggest the
following rule for learning rate $\eta_{k}$, when the denominator is positive:

$\eta_{k}^{2}\leq\frac{\left\|x_{k}-x_{k-1}\right\|^{2}}{3\left\|g_{k}\right\|^{2}-4\left\langle
g_{k},g_{k-1}\right\rangle}.$

Our method arose from attempts to adapt this learning rate to the stochastic
case.

##  3 Augmenting AdaGrad

For simplicity we consider the Euclidean step setting, with an unconstrained
domain. In this situation, as shown by Duchi et al. (2011), the AdaGrad update

$\displaystyle g_{k}$ $\displaystyle=\nabla f(x_{k},\xi_{k}),$ (5)
$\displaystyle x_{k+1}$ $\displaystyle=x_{k}-A_{k+1}^{-1}g_{k},$

obeys the inequality:

$\displaystyle
2\mathbb{E}\sum_{k=0}^{n}\left(f(x_{k})-f_{*}\right)\leq\left\|x_{0}-x_{*}\right\|_{A_{1}}^{2}+\mathbb{E}\sum_{k=1}^{n}\left\|x_{k}-x_{*}\right\|_{\left(A_{k+1}-A_{k}\right)}^{2}+\mathbb{E}\sum_{k=0}^{n}\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}.$
(6)

The right-hand side of this inequality consists of some initial conditions
(first term), an iterate distance error term (second term) and a gradient norm
term (third term). The key insight of the AdaGrad method is that the gradient
error term’s growth can be restricted to a $\mathcal{O}(\sqrt{n})$ rate over
time if the matrices $A_{k}$ are carefully chosen. For simplicity of notation,
consider the one-dimensional case for the remainder of this section. The
gradient error term grows each step by: $g_{k}^{2}/A_{k+1},$ and the overall
growth is controlled by using an inductive bound of the form:
$qA_{k}+g_{k}^{2}/A_{k+1}\leq qA_{k+1}.$ AdaGrad uses the choice
$A_{k+1}^{2}=\sum_{t=0}^{k}g_{t}^{2}$ and $q=2$ to keep the overall growth
rate $\mathcal{O}(\sqrt{n})$. In our approach, we add additional terms which
combine with the gradient error term. Inspired by hyper-gradient methods, we
use the following Lemma.

######  Lemma 1.

For any $A_{k}\succ 0$ and $\rho>0$, the sequences generated by (5) obey:

$0\leq-\mathbb{E}\left[\rho\left\langle
g_{k},g_{k-1}\right\rangle_{A_{k}^{-1}}|\mathcal{F}_{k}\right]+\rho\left[f(x_{k-1})-f(x_{k})\right].$

The first term is a scaling of the hyper-gradient; see Equation (4). Usually
this quantity is used to adapt the learning rate by applying stochastic
gradient descent on the learning rate, treating the hyper-parameter as a
parameter (hence the name). Our approach here is different, we introduce the
term into our convergence rate bound directly. The second term is a function
value difference between the current and previous points. This term is
particularly nice as it may be telescoped, and so it introduces an error term
that doesn’t grow over time. The degree of adaptivity is controlled by $\rho$,
which is a tunable parameter of our method. The theory we develop in Section
4.4 suggests that $\rho=2$ is the best default choice, although larger values
can lead to greater adaptivity.

Using Lemma 1 and the AdaGrad result (6), we may derive the following bound.

######  Theorem 2.

For $k\geq 1$, Algorithm 1 satisfies:

$\displaystyle
2\mathbb{E}\sum_{k=0}^{n}\left[f(x_{k})-f_{*}\right]\leavevmode\nobreak\
\leq\leavevmode\nobreak\ $
$\displaystyle\left\|x_{0}-x_{*}\right\|_{A_{1}}^{2}+\rho\left[f(x_{0})-f_{*}\right]\leavevmode\nobreak\
+\leavevmode\nobreak\
\mathbb{E}\sum_{k=1}^{n}\left\|x_{k}-x_{*}\right\|_{\left(A_{k+1}-A_{k}\right)}^{2}$
$\displaystyle+\mathbb{E}\sum_{k=0}^{n}\left[\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k-1},g_{k}\right\rangle_{A_{k}^{-1}}\right].$

Compared to the bound for AdaGrad in Equation (6), the gradient norm error
term now grows each step by a term:

$\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k},g_{k-1}\right\rangle_{A_{k}^{-1}}\quad\text{v.s.}\quad\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}.$

This error term is potentially much smaller, or even negative, when
consecutive gradients are highly positively correlated. This corresponds to
the situation in hyper-gradient methods where the learning rate should be
increased. This error term is larger in the opposite situation, where
consecutive gradients are negatively correlated, and the learning rate should
be decreased. In the “goldilocks” zone, the hyper-gradient term will wobble
around 0, and should not contribute significantly to the accumulated error
over time.

Ideally, we would like to apply a similar inductive bound to this error term
as was applied for AdaGrad:

$qA_{k}+\frac{\left\|g_{k}\right\|^{2}}{A_{k+1}}-\rho\frac{\left\langle
g_{k},g_{k-1}\right\rangle}{A_{k}}\leq qA_{k+1},$

with $A$ given by (3). This schema is the motivation for development of our
GradaGrad method. Applied directly we run into a number of problems that must
be solved to arrive at a practical method:

  1. i)

The above update can potentially yield the square root of a negative number.
The approach we develop in the following sections, shown in Figure 1, steadily
grows the learning rate, avoiding the divergence towards infinity encountered
in the naive implementation in Equation (3).

  2. ii)

When $A_{k}$ is allowed to both expand and shrink, generally, the iterate-
distance error term, i.e., second line in (6), may grow at a $\mathcal{O}(n)$
rate, rather than the $\mathcal{O}(\sqrt{n})$ rate required for a non-trivial
convergence rate bound. It’s this term that prevents the use of arbitrary
step-size sequences. Our approach maintains the $\mathcal{O}(\sqrt{n})$ rate
by construction.

  3. iii)

The hyper-gradient term is divided by $A_{k-1}$ not $A_{k}$, which
significantly breaks the bounding approach used by AdaGrad. The bound may be
violated both when the error term is positive and negative.

##  4 The GradaGrad Method

In this section, we show how to address the issues listed above through a
reparametrization technique, present the full version of GradaGrad with
momentum, and establish its convergence properties.

###  4.1 Controlling error through reparameterization

The key to controlling the error terms that arise in the GradaGrad method is
the use of reparameterization of the learning rate. Like AdaGrad, our learning
rate will take the form:

$A_{k+1}^{-1}=\frac{\gamma_{k+1}}{\sqrt{\alpha_{k+1}}},$

where $\alpha_{k+1}$ is updated from $\alpha_{k}$ each step by adding some
additional term $v_{k}$. In GradaGrad, we allow the numerator $\gamma_{k+1}$
to change over time, compared to the fixed numerator used in AdaGrad. The
purpose of this additional flexibility is to allow us to reparameterize the
learning rate before applying the AdaGrad like additive update to
$\alpha_{k}$. We still consider coordinate-wise updates of the form:

$\alpha_{k+1}=\alpha_{k}+v_{k},$ (7)

where $v_{k}=\left\|g_{k}\right\|^{2}$ in AdaGrad and
$v_{k}=\left\|g_{k}\right\|^{2}-\rho\left\langle g_{k},g_{k-1}\right\rangle$
in GradaGrad. However, we apply these updates after a reparameterization:

$\frac{\gamma_{k+1}}{\sqrt{\alpha_{k+1}^{\prime}}}=\frac{\gamma_{k}}{\sqrt{\alpha_{k}}},$
(8)

that leaves the learning rate the same but changes the effect of adding
$v_{k}$. The update then consists of
$\alpha_{k+1}=\alpha_{k+1}^{\prime}+v_{k}$. We choose our reparameterization
so that when $v_{k}$ is negative, $\alpha_{k+1}=\alpha_{k}$ regardless of the
value of $v_{k}$. Solving Equation (8) under this condition gives the update

$\displaystyle\gamma_{k+1}$
$\displaystyle=\gamma_{k}\sqrt{\frac{\alpha^{\prime}_{k+1}}{\alpha_{k}}}=\gamma_{k}\sqrt{\frac{\alpha_{k+1}-v_{k}}{\alpha_{k}}}=\gamma_{k}\sqrt{\frac{\alpha_{k}-v_{k}}{\alpha_{k}}}=\gamma_{k}\sqrt{1-\frac{v_{k}}{\alpha_{k}}}.$
(9)

This results in the behavior observed in Figure 1, where the learning rates
smoothly increases as $v_{k}$ becomes more negative. Note that this correction
is only used when $v_{k}$ is negative, for positive $v_{k}$ the standard
AdaGrad update in Equation (7) is used without any reparameterization. This
update is very well behaved, the learning rate matches the value from using
$\alpha_{k+1}=\alpha_{k}-v_{k}$ up to a first order approximation, as the
gradients match at $v_{k}=0$ (this can be seen in Figure 1). The ratio
$v_{k}/\alpha_{k}$ is typically very small at the later stages of
optimization, so its behavior differs from the naive variant only at the
beginning of optimization.

Input: $\gamma_{0}$, $\rho$, $x_{0}$, $\beta$, $G_{\infty}$, $D_{\infty}$

Initialize $z_{0}=x_{0}$

for $k=0$ to $n$ do

$g_{k}=\nabla f(x_{k},\xi_{k})$

for $i=1$ to $D$ do

if $k=0$ then

$v_{i,k}=G_{\infty}^{2}$

else if $\gamma_{i,k}=D_{\infty}$ then

$v_{i,k}=g_{i,k}^{2}$

else

$v_{i,k}=g_{i,k}^{2}-\rho g_{i,k}m_{i,k-1}$

end if

if $v_{i,k}\geq 0$ then

$\gamma_{i,k+1}=\gamma_{i,k}$

$\alpha_{i,k+1}=\alpha_{i,k}+v_{i,k}$

else

$r_{i,k}=\left(\frac{\rho m_{i,k-1}}{g_{i,k}}\right)^{2}-1$

$v_{i,k}=\max\left(v_{i,k},-r_{i,k}\alpha_{i,k}\right)$

$\gamma_{i,k+1}=\min\left(\gamma_{i,k}\sqrt{1-\frac{v_{i,k}}{\alpha_{i,k}}},\,D_{\infty}\right)$

$\alpha_{i,k+1}=\alpha_{i,k}$

end if

$A_{i,k+1}=\sqrt{\alpha_{i,k+1}}/\gamma_{i,k+1}$

end for

$z_{k+1}=\text{proj}_{\mathcal{X}}(z_{k}-A_{k+1}^{-1}g_{k})$

$x_{k+1}=\beta x_{k}+\left(1-\beta\right)z_{k+1}$

$m_{k}=A_{k+1}\left(x_{k}-x_{k+1}\right)$

end for

algocfGradaGrad

###  4.2 GradaGrad with momentum

The scalar version of GradaGrad is detailed in Algorithm 1. It’s straight-
forward to extend the method to per-coordinate adaptivity (diagonal scaling),
and to include the use of momentum and projection onto a convex set
$\mathcal{X}$. This full version is given in Algorithm 4.1. The full version
includes 2 other changes needed to facilitate the analysis: firstly, the
$v_{0}$ update uses a maximum element-wise gradient bound $G$ rather than the
observed gradient $g_{0}$. This kind of change is also needed by other AdaGrad
variants (such as AdaGrad-DA) for the theory and can usually be omitted from
practical implementations of the method.

Secondly, we derive an explicit expression for the hyper-parameter $\rho$
rather than using a fixed value. This is a trade-off, as it adds some
additional complexity to the method but reduces the number of hyper-
parameters. The necessity of this parameter is detailed in Section 4.3. We
establish the following $\mathcal{O}(1/\sqrt{n})$ convergence rate bound for
the diagonal scaling variant, the scalar variant can be analysed using the
same techniques.

######  Theorem 3.

Define $G_{\infty}=\max_{k}\left\|g_{k}\right\|_{\infty}$,
$D_{\infty}\geq\sup_{x\in\mathcal{X}}\left\|x-x_{*}\right\|_{\infty}$ then the
function value at the average iterate $\bar{x}$ of GradaGrad converges at a
$\mathcal{O}(1/\sqrt{n})$ rate:

$\displaystyle\mathbb{E}\left[f(\bar{x}_{n+1})-f_{*}\right]$
$\displaystyle\leq\frac{1}{2(n+1)}\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]+3\sqrt{2+\rho}\,\frac{dG_{\infty}D_{\infty}^{2}}{\gamma_{0}\sqrt{n+1}}.$

where $\rho>0$ is a hyper-parameter and $\beta\in[0,1)$ is the classical
momentum parameter.

###  4.3 The restricted increase constraint

As mentioned in point (iii), care must be taken to control the growth of the
gradient noise term:

$\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k-1},g_{k}\right\rangle_{A_{k}^{-1}}.$

The AdaGrad approach to bounding the accumulated error:

$qA_{k}+\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k-1},g_{k}\right\rangle_{A_{k}^{-1}}\leq qA_{k+1}.$

breaks down when $v_{k}=\left\|g_{k}\right\|^{2}-\rho\left\langle
g_{k},g_{k-1}\right\rangle$ is negative, as this inequality can become
impossible to satisfy for any choice of $A_{k}\prec A_{k-1}$, preventing us
from using increasing learning rate sequences. Instead, we will impose the
following restriction:

$\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k-1},g_{k}\right\rangle_{A_{k}^{-1}}\leq 0.$ (10)

This is not always satisfied, even when
$\left\|g_{k}\right\|^{2}-\rho\left\langle g_{k-1},g_{k}\right\rangle\leq 0$,
as the gradient norm square term may grow very large if the learning rate
$A_{k}^{-1}$ is allowed to increase significantly between steps. By
rearranging terms, we see that this bound is satisfied when the ratio of
$A_{k+1}$ to $A_{k}$ is controlled:

$A_{k+1}A_{k}^{-1}\succeq\frac{\left\|g_{k}\right\|^{2}}{\rho\left\langle
g_{k-1},g_{k}\right\rangle}.$

We satisfy this constraint by using the following update in the GradaGrad
method:

######  Lemma 4.

Inequality (10) is satisfied if we clip $v_{k}$ using:

$v_{k}=\max\left(v_{k},-r\alpha_{k}\right)\,\,\text{where}\,\,r=\left(\frac{\rho\left\langle
g_{k-1},g_{k}\right\rangle}{\left\|g_{k}\right\|^{2}}\right)^{2}-1.$

###  4.4 The $\rho$ hyper-parameter

There are two values of $\rho$ for which the update has special properties
worth further discussion. Consider the use of $\rho=1$. The $\alpha_{k}$
sequence (without reparameterization) can be written as

$\displaystyle\alpha_{n+1}$
$\displaystyle=g_{0}^{2}+\sum_{k=1}^{n}\left[g_{k}^{2}-g_{k}g_{k-1}\right]=\frac{1}{2}g_{0}^{2}+\frac{1}{2}g_{n}^{2}+\frac{1}{2}\sum_{k=0}^{n-1}\left(g_{k+1}-g_{k}\right)^{2}.$

Although individual steps may increase or decrease the rate, the overall
growth will be dominated by the summation term, and so GradaGrad with $\rho=1$
will not result in significant learning rate growth. The accumulation of
differences of gradients instead of gradients has been explored before, under
the name of optimistic online gradient descent (Rakhlin & Sridharan, 2013;
Bach & Levy, 2019; Antonakopoulos & Mertikopoulos, 2021).

This method is still of interest as it can be used as a drop-in replacement
for AdaGrad, without the need for the complexities of reparameterization that
larger $\rho$ require. It can behave quite differently than AdaGrad, for
instance on the absolute value function example in Figure 1, it will make much
more rapid progress, the learning rate initially stays constant rather than
decreasing as $g_{k+1}-g_{k}$ is zero until it overshoots the minimum, and
starts decreasing from there.

The case of $\rho=2$ is the most natural as it results in simplifications of
the constants in the bounds. We would recommend $\rho=2$ as the default choice
in the GradaGrad method, as it results in adaptivity without introducing
excessive instability during training.

Figure 2: GradaGrad matches the training accuracy of benchmark methods without
tuning (with default hyper-parameters $\gamma_{0}=1$ and $\rho=2$). Figure 3:
$\rho$ has a significant effect on the adaptability of Scalar GradaGrad on
large deep learning problems such as DLRM. Our recommended default value of
$\rho=2$ is highly stable. Larger values of $\rho$ give a higher degree of
adaptivity but may result in less stable training.

##  5 Experiments

We evaluated GradaGrad on 6 baseline problems widely used in the optimization
literature (Chang & Lin, 2011; Dua & Graff, 2017)): Glass, Letter, Protein,
USPS, Vehicle and Vowel. We tested a binary logistic regression model, with no
further data preprocessing, using the scaled versions of each dataset (when
applicable) retrieved from the LIBSVM dataset repository. We compare our
method against Adam, SGD and AdaGrad, where for each method the learning rate
is chosen via grid search on a power-of-2 grid to give the highest final
accuracy averaged over the last 10 epochs. This tuning is extremely
advantageous to these baselines methods, and provides a difficult benchmark to
beat. We plot the average of 10 runs with different random seeds, with a 2
standard-error range overlaid. As shown in Figure 2, GradaGrad is able to
match these baseline methods using the same default hyper-parameters of
$\gamma_{0}=1.0$, and $\rho=2$ across all problems.

To illustrate the effect of the $\rho$ hyperparameter, we performed
experiments on the DLRM (Naumov et al., 2019) model on the Criteo Kaggle
Display Advertising Challenge Dataset, using the scalar variant of GradaGrad.
This model has more than 300 million parameters and is representative of
models used in industry. Batch size 128 and embedding dimension 16 were used,
with no regularization. Whereas the diagonal scaling variant of GradaGrad
performs well with $\rho=2$, the scalar variant requires larger $\rho$ for
large or high noise problems. Figure 3 shows the results of using $\gamma_{0}$
equal to 1, 1e-2, 1e-4, 1e-6 with different $\rho$ values. Values of $\rho=2$
and $\rho=10$ show a significant lack of adaptivity, whereas for large $\rho$
each of these rates performs similarly.

##  6 Discussion

Although GradaGrad’s learning rate is more adaptive than that of AdaGrad, it
still inherits many of the same limitations. The “any-time” sequences that
have a square-root growth denominator give optimal rates up to a constant
factor, but they sometimes give noticeably worse rates in practice than more
heavily hand-tuned learning rate schedules such as MADGRAD. (Defazio &
Jelassi, 2021).

The GradaGrad method is also limited in its ability to adjust $\gamma$. The
updates only increase $\gamma$, so the initial $\gamma_{0}$ must be no larger
than the optimal $\gamma$. The GradaGrad update allows $\gamma$ to grow at an
exponential rate, which means that small initial values can potentially be
used. Care must be taken when initializing with $\gamma_{0}$ that are very
small as numerical issues may come into play. Values below roughly $10^{-6}$
will likely result in loss of precision and no adaptivity when using single
precision floating point.

## Conclusion

The GradaGrad method has promising theoretical and practical properties that
make it a viable drop-in replacement for the AdaGrad method. It exhibits
greater adaptivity to the gradient sequence as it is able to exploit
correlations between consecutive gradients, going beyond the simple gradient
magnitudes used in the AdaGrad method. There are many open problems relating
to the $v_{k}$ sequence we described. Can the same sequence be used in other
AdaGrad family methods? We believe that $v_{k}$ sequences of the form we
describe here could have widespread applications. Can the upper bound on
$\gamma$ be removed? It does not appear to be necessary in practice.

## References

  * Antonakopoulos & Mertikopoulos (2021) Antonakopoulos, K. and Mertikopoulos, P.  Adaptive first-order methods revisited: Convex minimization without lipschitz requirements.  In Beygelzimer, A., Dauphin, Y., Liang, P., and Vaughan, J. W. (eds.), _Advances in Neural Information Processing Systems_ , 2021.
  * Bach & Levy (2019) Bach, F. R. and Levy, K. Y.  A universal algorithm for variational inequalities adaptive to smoothness and noise.  In Beygelzimer, A. and Hsu, D. (eds.), _Conference on Learning Theory, COLT 2019_ , volume 99 of _Proceedings of Machine Learning Research_ , pp. 164–194. PMLR, 2019.
  * Baydin et al. (2018) Baydin, A. G., Cornish, R., Rubio, D. M., Schmidt, M., and Wood, F.  Online learning rate adaptation with hypergradient descent.  In _Proceedings of the Sixth International Conference on Learning Representations (ICLR)_ , Vancouver, Canada, 2018.
  * Chang & Lin (2011) Chang, C.-C. and Lin, C.-J.  Libsvm: A library for support vector machines.  _ACM transactions on intelligent systems and technology (TIST)_ , 2(3):27, 2011.
  * Defazio & Jelassi (2021) Defazio, A. and Jelassi, S.  Adaptivity without compromise: A momentumized, adaptive, dual averaged gradient method for stochastic optimization, 2021.
  * Delyon & Juditsky (1993) Delyon, B. and Juditsky, A.  Accelerated stochastic approximation.  _SIAM Journal on Optimization_ , 3(4):868–881, 1993.
  * Dua & Graff (2017) Dua, D. and Graff, C.  UCI machine learning repository, 2017.  URL http://archive.ics.uci.edu/ml.
  * Duchi et al. (2011) Duchi, J., Hazan, E., and Singer, Y.  Adaptive subgradient methods for online learning and stochastic optimization.  _Journal of Machine Learning Research_ , 12(61):2121–2159, 2011.
  * Jacobs (1988) Jacobs, R. A.  Increased rates of convergence through learning rate adaption.  _Neural Networks_ , 1:295–307, 1988.
  * Kesten (1958) Kesten, H.  Accelerated stochastic approximation.  _Annals of Mathematical Statistics_ , 29(1):41–59, 1958.
  * Kingma & Ba (2014) Kingma, D. P. and Ba, J.  Adam: A method for stochastic optimization.  _Proceedings of the 3rd International Conference on Learning Representations (ICLR)_ , 2014.
  * Li & Orabona (2019) Li, X. and Orabona, F.  On the convergence of stochastic gradient descent with adaptive stepsizes.  In _Proceedings of the 22nd International Conference on Artificial Intelligence and Statistics (AISTATS) 2019_ , 2019.
  * Loizou et al. (2021) Loizou, N., Vaswani, S., Laradji, I., and Lacoste-Julien, S.  Stochastic polyak step-size for sgd: An adaptive learning rate for fast convergence.  _Proceedings of the 24th International Conference on Artifi- cial Intelligence and Statistics (AISTATS) 2021_ , 2021.
  * Mahmood et al. (2012) Mahmood, A. R., Sutton, R. S., Degris, T., and Pilarski, P. M.  Tuning-free step-size adaption.  In _Proceedings of the IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_ , pp. 2121–2124, 2012.
  * Malitsky & Mishchenko (2019) Malitsky, Y. and Mishchenko, K.  Adaptive gradient descent without descent.  _arXiv preprint arXiv:1910.09529_ , 2019.
  * Mirzoakhmedov & Uryasev (1983) Mirzoakhmedov, F. and Uryasev, S. P.  Adaptive step adjustment for a stochastic optimization algorithm.  _Zh. Vychisl. Mat. Mat. Fiz._ , 23(6):1314–1325, 1983.  [U.S.S.R. Comput. Math. Math. Phys. 23:6, 1983].
  * Naumov et al. (2019) Naumov, M., Mudigere, D., Shi, H. M., Huang, J., Sundaraman, N., Park, J., Wang, X., Gupta, U., Wu, C., Azzolini, A. G., Dzhulgakov, D., Mallevich, A., Cherniavskii, I., Lu, Y., Krishnamoorthi, R., Yu, A., Kondratenko, V., Pereira, S., Chen, X., Chen, W., Rao, V., Jia, B., Xiong, L., and Smelyanskiy, M.  Deep learning recommendation model for personalization and recommendation systems.  _CoRR_ , abs/1906.00091, 2019.
  * Orabona & Tommasi (2017) Orabona, F. and Tommasi, T.  Training deep networks without learning rates through coin betting.  In _Advances in Neural Information Processing Systems_ , volume 30. Curran Associates, Inc., 2017.
  * Paquette & Scheinberg (2020) Paquette, C. and Scheinberg, K.  A stochastic line search method with expected coplexity analysis.  _SIAM Journal on Optimization_ , 30(1):349–376, 2020.
  * Rakhlin & Sridharan (2013) Rakhlin, A. and Sridharan, K.  Online learning with predictable sequences.  In _Advances in Neural Information Processing Systems_ , 2013.
  * Reddi et al. (2019) Reddi, S. J., Kale, S., and Kumar, S.  On the convergence of adam and beyond.  _e-Preprint arXiv:1904.09237_ , 2019.
  * Schaul et al. (2013) Schaul, T., Zhang, S., and LeCun, Y.  No more pesky learning rates.  In _Proceedings of the 30th International Conference on Machine Learning_ , volume 28 of _Proceedings of Machine Learning Research_ , pp. 343–351. PMLR, 2013.
  * Schraudolph (1999) Schraudolph, N. N.  Local gain adaptation in stochastic gradient descent.  In _Proceedings of Nineth International Conference on Artificial Neural Networks (ICANN)_ , pp. 569–574, 1999.
  * Sutton (1992) Sutton, R. S.  Adapting bias by gradient descent: An incremental version of Delta-Bar-Delta.  In _Proceedings of the Tenth National Conference on Artificial Intelligence (AAAI’92)_ , pp. 171–176. The MIT Press, 1992.
  * Tieleman & Hinton (2012) Tieleman, T. and Hinton, G.  Lecture 6.5-rmsprop: Divide the gradient by a running average of its recent magnitude.  _COURSERA: Neural networks for machine learning_ , 4(2):26–31, 2012.
  * van Erven & Koolen (2016) van Erven, T. and Koolen, W. M.  Metagrad: Multiple learning rates in online learning.  In _Advances in Neural Information Processing Systems_ , 2016.
  * Vaswani et al. (2019) Vaswani, S., Mishkin, A., Laradji, I., Schmidt, M., Gidel, G., and Lacoste-Julien, S.  Painless stochastic gradient: Interpolation, line-search, and convergence rates.  In _Advances in Neural Information Processing Systems_ , volume 32. Curran Associates, Inc., 2019.
  * Zhang et al. (2020) Zhang, P., Lang, H., Liu, Q., and Xiao, L.  Statistical adaptive stochastic gradient methods.  e-Preprint: arXiv:2002.10597, 2020.

##  Appendix A Convergence Theory

We define the ancillary quantity $z^{\prime}$, so that the update in Algorithm
4.1 may be written as

$\displaystyle z_{k+1}^{\prime}$ $\displaystyle=z_{k}-A_{k+1}^{-1}g_{k},$
$\displaystyle z_{k+1}$
$\displaystyle=\text{proj}_{\mathcal{X}}(z_{k+1}^{\prime}),$ $\displaystyle
x_{k+1}$ $\displaystyle=\beta x_{k}+\left(1-\beta\right)z_{k+1}.$ (11)

Then define $m_{k}$ via

$x_{k+1}=x_{k}-A_{k+1}^{-1}m_{k}.$

Note that $m_{k}=A_{k}\left(x_{k}-x_{k+1}\right)$. Note also that:

$z_{k}=\frac{1}{1-\beta}x_{k}-\frac{\beta}{1-\beta}x_{k-1}.$

######  Lemma 5.

For $k\geq 1$:

$\mathbb{E}_{\xi_{k}}\left\langle
g_{k},m_{k-1}\right\rangle_{A_{k}^{-1}}\leq\left[f(x_{k-1})-f(x_{k})\right].$

###### Proof.

Using $x_{k}=x_{k-1}-A_{k}^{-1}m_{k-1}$ ,

$\displaystyle\mathbb{E}_{\xi_{k}}\left\langle
g_{k},m_{k-1}\right\rangle_{A_{k}^{-1}}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\langle
g_{k},x_{k-1}-x_{k}\right\rangle$ $\displaystyle=\left\langle\nabla
f(x_{k}),x_{k-1}-x_{k}\right\rangle$ $\displaystyle\leq f(x_{k-1})-f(x_{k}).$

∎

######  Lemma 6.

Consider $k\geq 1$ and $\rho\geq 1$. Then the iterate sequence generated by
Equation (11) obeys the following bound:

$\displaystyle 2\left[f(x_{k})-f_{*}\right]$
$\displaystyle\leq\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-\mathbb{E}_{\xi_{k}}\left\|z_{k+1}-x_{*}\right\|_{A_{k+1}}^{2}$
$\displaystyle-\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{k})-f_{*}\right]+\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{k-1})-f_{*}\right]$
$\displaystyle+\mathbb{E}_{\xi_{k}}\left[\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k},m_{k-1}\right\rangle_{A_{k}^{-1}}\right].$

###### Proof.

We start by expanding $x_{k+1}$ around $x_{*}$.

$\displaystyle\mathbb{E}_{\xi_{k}}\left\|z_{k+1}-x_{*}\right\|_{A_{k+1}}^{2}$
$\displaystyle\leq\mathbb{E}_{\xi_{k}}\left\|z_{k+1}^{\prime}-x_{*}\right\|_{A_{k+1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}+z_{k+1}^{\prime}-z_{k}\right\|_{A_{k+1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}+2\mathbb{E}_{\xi_{k}}\left\langle
A_{k+1}\left(z_{k+1}^{\prime}-z_{k}\right),z_{k}-x_{*}\right\rangle+\mathbb{E}_{\xi_{k}}\left\|z_{k+1}^{\prime}-z_{k}\right\|_{A_{k+1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-2\mathbb{E}_{\xi_{k}}\left\langle
g_{k},z_{k}-x_{*}\right\rangle+\mathbb{E}\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-2\left\langle\nabla
f(x_{k}),\frac{1}{1-\beta}x_{k}-\frac{\beta}{1-\beta}x_{k-1}-x_{*}\right\rangle+\mathbb{E}_{\xi_{k}}\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-2\left\langle\nabla
f(x_{k}),\frac{\beta}{1-\beta}x_{k}-\frac{\beta}{1-\beta}x_{k-1}+x_{k}-x_{*}\right\rangle+\mathbb{E}_{\xi_{k}}\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-2\left\langle\nabla
f(x_{k}),x_{k}-x_{*}\right\rangle-\frac{2\beta}{1-\beta}\left\langle\nabla
f(x_{k}),x_{k}-x_{k-1}\right\rangle+\mathbb{E}_{\xi_{k}}\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}$
$\displaystyle\leq\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-2\left[f(x_{k})-f_{*}\right]-\frac{2\beta}{1-\beta}\left[f(x_{k})-f(x_{k-1})\right]+\mathbb{E}_{\xi_{k}}\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}$

We combine with Lemma 5 multiplied by $\rho$:

$\displaystyle\mathbb{E}_{\xi_{k}}\left\|z_{k+1}-x_{*}\right\|_{A_{k+1}}^{2}$
$\displaystyle\leq\mathbb{E}_{\xi_{k}}\left\|z_{k}-x_{*}\right\|_{A_{k+1}}^{2}-2\left[f(x_{k})-f_{*}\right]+\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{k-1})-f(x_{k})\right]$
$\displaystyle+\mathbb{E}_{\xi_{k}}\left[\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k},m_{k-1}\right\rangle_{A_{k}^{-1}}\right].$

∎

###  A.1 Proof of Theorem 2

We prove a more general form of Theorem 2, covering the use of momentum and
projection. The variant listed in the body of the paper is the special case
where $\mathcal{X}=\mathbb{R}^{d}$ and $\beta=0$. Consider the case where
$k=0$. Then:

$\displaystyle\mathbb{E}_{\xi_{0}}\left\|z_{1}-x_{*}\right\|_{A_{1}}^{2}$
$\displaystyle\leq\mathbb{E}_{\xi_{0}}\left\|z_{1}^{\prime}-x_{*}\right\|_{A_{1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{0}}\left\|z_{0}-x_{*}+z_{1}^{\prime}-z_{0}\right\|_{A_{1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{0}}\left\|z_{0}-x_{*}\right\|_{A_{1}}^{2}+2\left\langle
A_{1}\left(z_{1}^{\prime}-z_{0}\right),z_{0}-x_{*}\right\rangle+\mathbb{E}_{\xi_{0}}\left\|z_{1}^{\prime}-z_{0}\right\|_{A_{1}}^{2}$
$\displaystyle=\mathbb{E}_{\xi_{0}}\left\|z_{0}-x_{*}\right\|_{A_{1}}^{2}-2\mathbb{E}\left\langle
g_{0},z_{0}-x_{*}\right\rangle+\mathbb{E}_{\xi_{0}}\left\|g_{0}\right\|_{A_{1}^{-1}}^{2}$
$\displaystyle\leq\mathbb{E}_{\xi_{0}}\left\|z_{0}-x_{*}\right\|_{A_{1}}^{2}-2\left[f(x_{0})-f_{*}\right]+\mathbb{E}_{\xi_{0}}\left\|g_{0}\right\|_{A_{1}^{-1}}^{2}.$

Using the law of total expectation, we telescope Lemma 6 together with the
base case from $k=0$ to $n$:

$\displaystyle 2\mathbb{E}\sum_{k=0}^{n}\left[f(x_{k})-f_{*}\right]$
$\displaystyle\leq\mathbb{E}\left\|z_{0}-x_{*}\right\|_{A_{1}}^{2}+\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+\mathbb{E}\sum_{k=1}^{n}\left\|z_{k}-x_{*}\right\|_{\left(A_{k+1}-A_{k}\right)}^{2}$
$\displaystyle+\mathbb{E}\sum_{k=0}^{n}\left[\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k},m_{k-1}\right\rangle_{A_{k}^{-1}}\right].$

######  Lemma 7.

Define $D_{\infty}=\sup_{x\in\mathcal{X}}\left\|x-x_{*}\right\|$. The
GradaGrad update satisfies:

$\left\|z_{0}-x_{*}\right\|_{A_{1}}^{2}+\mathbb{E}\sum_{k=1}^{n}\left[\left\|z_{k}-x_{*}\right\|_{\left(A_{k+1}-A_{k}\right)}^{2}\right]\leq\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sum_{i=1}^{d}\sqrt{\alpha_{i,n+1}}.$

###### Proof.

We consider the 1D case for simplicity as the bound is separable in the
dimension of the problem. This allows us to simplify the left hand side to

$\left\|z_{0}-x_{*}\right\|_{A_{0}}^{2}+\mathbb{E}\sum_{k=1}^{n}\left[\left(\frac{\sqrt{\alpha_{k+1}}}{\gamma_{k+1}}-\frac{\sqrt{\alpha_{k}}}{\gamma_{k}}\right)\left\|z_{k}-x_{*}\right\|^{2}\right].$

We will prove this by induction. Consider the base case, $v_{0}$ is always
positive so:

$\left\|z_{0}-x_{*}\right\|_{A_{0}}^{2}\leq\frac{D_{\infty}^{2}}{\gamma_{0}}\sqrt{v_{0}}\leq\frac{D_{\infty}^{2}}{\gamma_{0}}\sqrt{\alpha_{1}}.$

Now consider $k>0$. We have two cases corresponding to the sign of $v_{k}$. If
$v_{k}$ is negative, then $\gamma_{k+1}$ and $\alpha_{k+1}$ are rescaled, and
$A_{k}$ decreases from $A_{k-1}$ so
$\frac{\sqrt{\alpha_{k+1}}}{\gamma_{k+1}}-\frac{\sqrt{\alpha_{k}}}{\gamma_{k}}$
is negative, and so:

$\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sqrt{\alpha_{k}}+\mathbb{E}\left(\frac{\sqrt{\alpha_{k+1}}}{\gamma_{k+1}}-\frac{\sqrt{\alpha_{k}}}{\gamma_{k}}\right)\left\|z_{k}-x_{*}\right\|^{2}\leq\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sqrt{\alpha_{k+1}}.$

Now consider the remaining case, where $v_{k}$ is positive, in which case
$\gamma$ increases and $\alpha$ stays unchanged:

$\displaystyle\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sqrt{\alpha_{k}}+\mathbb{E}\left(\frac{\sqrt{\alpha_{k+1}}}{\gamma_{k+1}}-\frac{\sqrt{\alpha_{k}}}{\gamma_{k}}\right)\left\|z_{k}-x_{*}\right\|^{2}$
$\displaystyle\leq\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sqrt{\alpha_{k}}$
$\displaystyle=\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sqrt{\alpha_{k+1}}$

∎

###  A.2 Gradient error bounding

######  Lemma 8.

When $v_{k}\leq 0$:

$\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
m_{k-1},g_{k}\right\rangle_{A_{k}^{-1}}\leq 0,$

when

$v_{k}\geq-r\alpha_{k},$ $r=\frac{1}{h_{k}^{2}}-1,$
$h_{k}=\frac{g_{k}^{2}}{\rho m_{k-1}g_{k}}.$

We again consider the 1D case as the inequality is separable in the dimension.
So the required inequality is:

$\frac{g_{k}^{2}}{A_{k+1}}-\rho\frac{m_{k-1}g_{k}}{A_{k}}\leq 0.$ (12)

Rearranging:

$g_{k}^{2}\leq\rho\frac{A_{k+1}}{A_{k}}m_{k-1}g_{k},$
$\frac{A_{k+1}}{A_{k}}\geq\frac{g_{k}^{2}}{\rho m_{k-1}g_{k}},$

Define:

$h_{k}=\frac{g_{k}^{2}}{\rho m_{k-1}g_{k}}\geq 0.$

We want to find a value of $r$ such that:

$\frac{A_{k+1}}{A_{k}}\geq h_{k}.$

Using $v_{k}=-r_{k}\alpha_{k}$ we have:

$\frac{\frac{\sqrt{\alpha_{k+1}}}{\gamma_{k+1}}}{\frac{\sqrt{\alpha_{k}}}{\gamma_{k}}}\geq
h_{k},$
$\frac{\frac{\sqrt{\alpha_{k}}}{\gamma_{k}\sqrt{\frac{\alpha_{k}-v_{k}}{\alpha_{k}}}}}{\frac{\sqrt{\alpha_{k}}}{\gamma_{k}}}\geq
h_{k},$

$\frac{1}{\sqrt{\frac{\alpha_{k}-v_{k}}{\alpha_{k}}}}\geq h_{k},$
$\frac{1}{h_{k}}\geq\sqrt{\left(1+r\right)},$ $r\leq\frac{1}{h_{k}^{2}}-1.$

Therefore any $r$ value smaller than $1/h_{k}^{2}-1$ is sufficient.

######  Lemma 9.

Consider the 1D case, with $k\geq 1$. Let $q\geq
2\sqrt{2+\rho}\gamma_{k+1}^{2}$ if $v_{k}\geq 0$:

$qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho g_{k}m_{k-1}}{A_{k}}\leq
qA_{k+1}.$

###### Proof.

We start by applying a case argument. Consider when $-g_{k}m_{k-1}\geq 0$.
From the fact that $\alpha_{0}=G_{\infty}^{2}$ we have
$\frac{1}{A_{k}}\leq\frac{\sqrt{2+\rho}}{A_{k+1}}$, and therefore:

$\displaystyle qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k}}$ $\displaystyle\leq
qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\sqrt{2+\rho}\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}$ $\displaystyle\leq
qA_{k}+\sqrt{2+\rho}\frac{g_{k}^{2}}{A_{k+1}}-\sqrt{2+\rho}\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}.$

In the other case, where $-g_{k}m_{k-1}<0$, then we can instead use
$A_{k+1}\geq A_{k}$ and thus $-\frac{1}{A_{k}}\leq-\frac{1}{A_{k+1}}$ (given
that $v_{k}$ is assumed positive), so:

$\displaystyle qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k}}$ $\displaystyle\leq
qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho g_{k}m_{k-1}}{A_{k+1}}$
$\displaystyle\leq
qA_{k}+\sqrt{2+\rho}\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}\right].$

Using the notation $A_{k}=\gamma_{k+1}^{-1}\sqrt{b_{k}-\alpha_{k}}$, where
$b_{k}=\sum_{t=0}^{k}\left(g_{t}^{2}-\rho g_{t}m_{t-1}\right)$ and
$\alpha_{k}=g_{k}^{2}-\rho g_{k}m_{k-1}$. We have from the concavity of the
square-root function that:

$\displaystyle\sqrt{b_{k}-\alpha_{k}}$
$\displaystyle\leq\sqrt{b_{k}}+\frac{b_{k}-\alpha_{k}-b_{k}}{2\sqrt{b_{k}}}$
$\displaystyle=\sqrt{b_{k}}-\frac{\alpha_{k}}{2\sqrt{b_{k}}}.$

Therefore:

$\displaystyle qA_{k}+\sqrt{2+\rho}\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}\right]$ $\displaystyle\leq
qA_{k+1}-\frac{q}{2\gamma_{k+1}}\left[\frac{g_{k}^{2}}{\sqrt{b_{k}}}-\frac{\rho
g_{k}m_{k-1}}{\sqrt{b_{k}}}\right]+\sqrt{2+\rho}\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}\right]$
$\displaystyle=qA_{k+1}-\frac{q}{2\gamma_{k+1}^{2}}\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}\right]+\sqrt{2+\rho}\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}\right]$
$\displaystyle=qA_{k+1}+\left(-\frac{2\sqrt{2+\rho}\gamma_{k+1}^{2}}{2\gamma_{k+1}^{2}}+\sqrt{2+\rho}\right)\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k+1}}\right]$ $\displaystyle=qA_{k+1}.$

∎

######  Lemma 10.

Let $q=2\sqrt{2+\rho}\gamma_{n}^{2}.$ Then:

$\sum_{k=0}^{n}\left[\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k-1},g_{k}\right\rangle_{A_{k}^{-1}}\right]\leq
q\sum_{i=1}^{d}A_{i,1}+\frac{q}{\gamma_{0}}\sum_{i=1}^{d}\sqrt{\alpha_{i,n+1}}.$

###### Proof.

Without loss of generality, we again consider the 1D case. Recall that if
$v_{k}<0$ (which can only occur when $k\geq 1$):

$\displaystyle qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k}}$ $\displaystyle\leq qA_{k}$
$\displaystyle=qA_{k+1}+q\left(A_{k}-A_{k+1}\right),$

and if $v_{k}>0$ and $k\geq 1$:

$qA_{k}+\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho g_{k}m_{k-1}}{A_{k}}\leq
qA_{k+1}.$

We will telescope this bound from the base-case $\frac{g_{0}^{2}}{A_{1}}\leq
qA_{1}$ (which is trivially true) to $A_{n+1}$. Combining we have:

$\displaystyle\sum_{k=0}^{n}\left[\frac{g_{k}^{2}}{A_{k+1}}-\frac{\rho
g_{k}m_{k-1}}{A_{k}}\right]$ $\displaystyle\leq
qA_{n+1}+q\sum_{k=1}^{n}I\left[v_{k}\leq 0\right]\left(A_{k}-A_{k+1}\right)$
$\displaystyle=qA_{1}-q\sum_{k=1}^{n}\left(A_{k}-A_{k+1}\right)+q\sum_{k=1}^{n}I\left[v_{k}\leq
0\right]\left(A_{k}-A_{k+1}\right)$
$\displaystyle=qA_{1}-q\sum_{k=1}^{n}I\left[v_{k}\geq
0\right]\left(A_{k}-A_{k+1}\right)$
$\displaystyle=qA_{1}+q\sum_{k=1}^{n}I\left[v_{k}>0\right]\left(A_{k+1}-A_{k}\right)$
$\displaystyle=qA_{1}+q\sum_{k=1}^{n}I\left[v_{k}>0\right]\frac{1}{\gamma_{k}}\left(\sqrt{\alpha_{k+1}}-\sqrt{\alpha_{k}}\right)$
$\displaystyle\leq
qA_{1}+\frac{q}{\gamma_{0}}\sum_{k=1}^{n}I\left[v_{k}>0\right]\left(\sqrt{\alpha_{k+1}}-\sqrt{\alpha_{k}}\right)$
$\displaystyle=qA_{1}+\frac{q\sqrt{\alpha_{n+1}}}{\gamma_{0}}.$

∎

###  A.3 Proof of Theorem 3

Recall the general form of Theorem 2 from Section A.1:

$\displaystyle 2\mathbb{E}\sum_{k=0}^{n}\left[f(x_{k})-f_{*}\right]$
$\displaystyle\leq\mathbb{E}\left\|z_{0}-x_{*}\right\|_{A_{1}}^{2}+\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+\mathbb{E}\sum_{k=1}^{n}\left\|z_{k}-x_{*}\right\|_{\left(A_{k+1}-A_{k}\right)}^{2}$
$\displaystyle+\mathbb{E}\sum_{k=0}^{n}\left[\left\|g_{k}\right\|_{A_{k+1}^{-1}}^{2}-\rho\left\langle
g_{k},m_{k-1}\right\rangle_{A_{k}^{-1}}\right].$

We apply the iterate error bound (Lemma 7) and the gradient error bound (Lemma
10) to give:

$\displaystyle 2\mathbb{E}\sum_{k=0}^{n}\left[f(x_{k})-f_{*}\right]$
$\displaystyle\leq\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+\frac{D_{\infty}^{2}}{\gamma_{0}}\mathbb{E}\sum_{i=1}^{d}\sqrt{\alpha_{i,n+1}}+q\mathbb{E}\sum_{i=1}^{d}A_{i,1}+\mathbb{E}\frac{q}{\gamma_{0}}\sum_{i=1}^{d}\sqrt{\alpha_{i,n+1}}.$

Now using
$f(\bar{x}_{n+1})-f_{*}\leq\frac{1}{n+1}\sum_{k=0}^{n}\left[f(x_{k})-f_{*}\right]$,
we have:

$\displaystyle\mathbb{E}\left[f(\bar{x}_{n+1})-f_{*}\right]$
$\displaystyle\leq\frac{1}{2(n+1)}\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+\left(\frac{D_{\infty}^{2}}{2\gamma_{0}(n+1)}+\frac{q}{2\gamma_{0}(n+1)}\right)\mathbb{E}\sum_{i=1}^{d}\sqrt{\alpha_{i,n+1}}+\frac{q}{2(n+1)}\mathbb{E}\sum_{i=1}^{d}A_{i,1}.$

Now using that $\alpha_{i,n+1}\leq(n+1)G_{\infty}^{2},$ and $A_{i1}\leq
G_{\infty}/\gamma_{0}$ and the definition of $q$:

$\displaystyle\mathbb{E}\left[f(\bar{x}_{n+1})-f_{*}\right]$
$\displaystyle\leq\frac{1}{2(n+1)}\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+\left(\frac{D_{\infty}^{2}}{2\gamma_{0}(n+1)}+\frac{q}{2\gamma_{0}(n+1)}\right)d\sqrt{n+1}G_{\infty}+\frac{qd}{2(n+1)\gamma_{0}}G_{\infty}.$

Simplifying:

$\displaystyle\mathbb{E}\left[f(\bar{x}_{n+1})-f_{*}\right]$
$\displaystyle\leq\frac{1}{2(n+1)}\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+\left(\frac{D_{\infty}^{2}}{2}+q\right)\frac{dG_{\infty}}{\gamma_{0}\sqrt{n+1}}.$

Now further using that $q=2\sqrt{2+\rho}\gamma_{n}^{2}\leq
2\sqrt{2+\rho}D_{\infty}^{2}$:

$\displaystyle\mathbb{E}\left[f(\bar{x}_{n+1})-f_{*}\right]$
$\displaystyle\leq\frac{1}{2(n+1)}\left(\frac{2\beta}{1-\beta}+\rho\right)\left[f(x_{0})-f_{*}\right]$
$\displaystyle+3\sqrt{2+\rho}\frac{dG_{\infty}D_{\infty}^{2}}{\gamma_{0}\sqrt{n+1}}.$

##  Appendix B Further experimental details

For the convex experiments, we trained in PyTorch using CPU training. Each
training run took several minutes at most, depending on the dataset. Our DLRM
training was run on a V100 GPU, and took approximately 6 hours per training
run.