<figure><img src="image_1.png" /></figure>

HAL
open science

# Tracking solutions of time-varying variational inequalities

Hédi Hadiji, Sarah Sachs, Cristóbal Guzmán

► To cite this version:

Hédi Hadiji, Sarah Sachs, Cristóbal Guzmán. Tracking solutions of time-varying variational inequalities. 2024. hal-04617351

HAL Id: hal-04617351

https://hal.science/hal-04617351

Preprint submitted on 19 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# Tracking solutions of time-varying variational inequalities

Hédi Hadiji, Sarah Sachs, Cristóbal Guzmán

► To cite this version:

Hédi Hadiji, Sarah Sachs, Cristóbal Guzmán. Tracking solutions of time-varying variational inequalities. 2025. hal-04617351v2

HAL Id: hal-04617351

https://hal.science/hal-04617351v2

Preprint submitted on 3 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# TRACKING SOLUTIONS OF TIME-VARYING VARIATIONAL INEQUALITIES

Hedi Hadiji*¹, Sarah Sachs*² and Cristóbal Guzmán²

¹Laboratoire des signaux et systèmes, Univ. Paris-Saclay, CNRS, CentraleSupélec

²Korteweg-de Vries Institute for Mathematics, University of Amsterdam

³Pontificia Universidad Católica de Chile, Institute for Mathematical and Computational Engineering, Facultad de Matemáticas and Escuela de Ingeniería

*equal contribution

June 19, 2024

## ABSTRACT

Tracking the solution of time-varying variational inequalities is an important problem with applications in game theory, optimization, and machine learning. Existing work considers time-varying games or time-varying optimization problems. For strongly convex optimization problems or strongly monotone games, these results provide tracking guarantees under the assumption that the variation of the time-varying problem is restrained, that is, problems with a sublinear solution path. In this work we extend existing results in two ways: In our first result, we provide tracking bounds for (1) variational inequalities with a sublinear solution path but not necessarily monotone functions, and (2) for periodic time-varying variational inequalities that do not necessarily have a sublinear solution path-length. Our second main contribution is an extensive study of the convergence behavior and trajectory of discrete dynamical systems of periodic time-varying VI. We show that these systems can exhibit provably chaotic behavior or can converge to the solution. Finally, we illustrate our theoretical results with experiments.

## 1 Introduction and setting.

Solutions of variational inequalities generalize minimizers of functions, saddle-points, fixed points, and many notions of equilibria across various disciplines. Practitioners often encounter problems in which the parameters of the variational inequality have a component that depends on time. We study the problem of tracking the sequence of solutions in time-varying variational inequalities: that is, sequentially selecting points that are as close to the solution as possible.

Let $\mathcal{Z}$ be a closed convex subset of $\mathbb{R}^d$ for $d \ge 1$. We let $\mathcal{C}(\mathcal{Z}, \mathbb{R}^d)$ denote the set of all continuous functions $\mathcal{Z} \to \mathbb{R}^d$. The finite-dimensional Variational Inequality Problem (VI problem) for the operator $F$ on the domain $\mathcal{Z}$, referred to as VIP($F, \mathcal{Z}$) from now on, is defined as: Given $F \in \mathcal{C}(\mathcal{Z}, \mathbb{R}^d)$,

$$ \text{find } Z^* \in \mathcal{Z} \text{ such that } \langle F(Z^*), Z - Z^* \rangle \ge 0 \quad \text{for all } Z \in \mathcal{Z}. \qquad (\text{VIP}) $$

For $\mu > 0$, the operator $F$ is said to be $\mu$-strongly monotone if

$$ \forall Z, Z' \in \mathcal{Z}: \quad \langle F(Z) - F(Z'), Z - Z' \rangle \ge \mu \|Z - Z'\|^2. \qquad (\mu\text{-Mono}) $$

The operator $F$ is monotone if the above holds for $\mu = 0$. Given $L \ge 0$, we say $F$ is $L$-Lipschitz if

$$ \forall Z, Z' \in \mathcal{Z}: \quad \|F(Z) - F(Z')\| \le L \|Z - Z'\| \qquad (L\text{-Lip}) $$

Many common problems in mathematics can be formulated as VI problems. Convex optimization is a prime example: Let $\mathcal{Z}$ be a closed convex set. Given a differentiable and convex function $f: \mathcal{Z} \to \mathbb{R}$, if we set $F = \nabla f$ then $u^*$ is a solution of VIP($F, \mathcal{Z}$) if and only if

$$ u^* \in \operatorname*{argmin}_{x \in \mathcal{Z}} f(x). \qquad (\text{OPT}) $$
# TRACKING SOLUTIONS OF TIME-VARYING VARIATIONAL INEQUALITIES

Hedi Hadiji*¹, Sarah Sachs*² and Cristóbal Guzmán²

¹Laboratoire des signaux et systèmes, Univ. Paris-Saclay, CNRS, CentraleSupélec

²Korteweg-de Vries Institute for Mathematics, University of Amsterdam

³Pontificia Universidad Católica de Chile, Institute for Mathematical and Computational Engineering, Facultad de Matemáticas and Escuela de Ingeniería

*equal contribution

November 3, 2025

## ABSTRACT

Tracking the solution of time-varying variational inequalities is an important problem with applications in game theory, optimization, and machine learning. Existing work considers time-varying games or time-varying optimization problems. For strongly convex optimization problems or strongly monotone games, such results provide tracking guarantees under the assumption that the variation of the time-varying problem is restrained, that is, problems with a sublinear solution path. We extend existing results in two ways: In our first result, we provide tracking bounds for (1) variational inequalities with a sublinear solution path but not necessarily monotone functions, and (2) for periodic time-varying variational inequalities that do not necessarily have a sublinear solution path length. Our second main contribution is an extensive study of the convergence behavior and trajectory of discrete dynamical systems of periodic time-varying VI. We show that these systems can either exhibit provably chaotic behavior or can converge to the solution.

## 1 Introduction and setting.

Solutions of variational inequalities generalize minimizers of functions, saddle-points, fixed points, and many notions of equilibria across various disciplines. Practitioners often encounter problems in which the parameters of the variational inequality have a component that depends on time. We study the problem of tracking the sequence of solutions in time-varying variational inequalities in an online setting: that is, sequentially predicting the solutions of a sequence of time-varying variational inequality problems.

**Formal Setting** Let $\mathcal{Z}$ be a closed convex subset of $\mathbb{R}^d$ for $d \ge 1$. We let $\mathcal{C}(\mathcal{Z}, \mathbb{R}^d)$ denote the set of all continuous functions $\mathcal{Z} \to \mathbb{R}^d$. The finite-dimensional Variational Inequality Problem (VI problem) for the operator $F$ on the domain $\mathcal{Z}$, referred to as VIP($F, \mathcal{Z}$) from now on, is defined as: Given $F \in \mathcal{C}(\mathcal{Z}, \mathbb{R}^d)$,

$$ \text{find } Z^* \in \mathcal{Z} \text{ such that } \langle F(Z^*), Z - Z^* \rangle \ge 0 \quad \text{for all } Z \in \mathcal{Z}. \qquad (\text{VIP}) $$

Let $\| \cdot \|$ denote any norm in $\mathbb{R}^d$. For $\mu > 0$, the operator $F$ is said to be $\mu$-strongly monotone if

$$ \forall Z, Z' \in \mathcal{Z}: \quad \langle F(Z) - F(Z'), Z - Z' \rangle \ge \mu \|Z - Z'\|^2. \qquad (\mu\text{-Mono}) $$

The operator $F$ is monotone if the above holds for $\mu = 0$. Given $L \ge 0$, we say $F$ is $L$-Lipschitz if

$$ \forall Z, Z' \in \mathcal{Z}: \quad \|F(Z) - F(Z')\| \le L \|Z - Z'\| \qquad (L\text{-Lip}) $$

Many common problems in mathematics can be formulated as VI problems. Prime examples are convex optimization and Nash equilibrium computation in concave games. See e.g., Rockafellar and Wets [1998, Chapter 12], Kinderlehrer
If the function f is μ-strongly convex and L-smooth, then the corresponding VI problem is μ-strongly monotone and L-Lipschitz. However, VI problems are far more general than optimization problems. To illustrate the relation to game theory, consider an n-player normal-form game with convex loss functions: Let $Z^{(i)}$ denote the (convex) action set of player i ∈ [n], and define $Z = \prod_{i \in [n]} Z^{(i)}$, and $Z^{(-i)} = \prod_{i \in [n] \setminus \{i\}} Z^{(i)}$. Finally, denote by $\nu^{(i)}: Z^{(i)} \times Z^{(-i)} \to \mathbb{R}$ the i-th player's loss function and assume that for all $x^{(-i)} \in Z^{(-i)}$, the loss $\nu^{(i)}(\cdot, x^{(-i)})$ is convex and differentiable. A point $\bar{x} = [\bar{x}^{(1)}, \dots, \bar{x}^{(n)}] \in Z$ is a Nash equilibrium if

$$ \forall i \in [n], \forall x^{(i)} \in Z^{(i)}: \quad \nu^{(i)}(x^{(i)}, \bar{x}^{(-i)}) \ge \nu^{(i)}(\bar{x}^{(i)}, \bar{x}^{(-i)}). \qquad (NE) $$

Consider the operator F mapping $(x^{(1)}, \dots, x^{(n)})$ to $[\nabla_{x^{(1)}} \nu^{(1)}(x^{(1)}, x^{(-1)}), \dots, \nabla_{x^{(n)}} \nu^{(n)}(x^{(n)}, x^{(-n)})]$, also called the pseudo-gradient of the game $(\nu^{(i)})$. The set of Nash equilibria is exactly the set of solutions of VIP($F$, $Z$). See, e.g., Rockafellar and Wets [1998, Chapter 12], or Kinderlehrer and Stampacchia [2000, Chapter 7], for more details and examples.

We study the problem of online tracking of solutions for time-varying VI. More precisely, we consider a sequence of operators $(F_t)_{t \in \mathbb{N}}$ revealed in an online fashion to a learner. At time step t, given the past observations, the learner outputs a candidate solution $Z_t$, and then observes the operator $F_t$. Formally, an online algorithm is a sequence of update rules $A_t: Z \times \mathcal{F}^{t-1} \to Z$ which operate online: at time step t the output of the algorithm may depend on $F_1, \dots, F_{t-1}, Z_{t-1}$, but not on $(F_s)_{s \ge t}$. An instance of a time-varying VI problem is characterized by the sequence $(F_t)$ together with the domain $Z$; we denote it by VIP$((F_t)_{t \in [T]}, Z)$.

We focus especially on sequences of operators that all admit a unique solution, denoted by $(Z_t^*)$. A sufficient, but not necessary, condition for uniqueness is that the operator be strongly monotone. For such problems, we examine the tracking guarantees of algorithms. Denoting the t-th output of algorithm A by $Z_t = A_t(Z_{t-1}, (F_s)_{s \in [t-1]})$, the tracking error of A is

$$ \tau_T(A) = \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2. $$

We study conditions under which this tracking error is sublinear, i.e., when there exists an $\alpha \in [0, 1)$ and constant $c \ge 0$, such that $\tau_T(A) \le cT^\alpha$. Note that guaranteeing a sublinear tracking error for time-varying VI problems may be impossible for some sequences of problems, especially for problems where the solutions vary too much. Indeed, there are simple instances such that for any online algorithm A and $T \in \mathbb{N}$, the tracking error grows linearly, as $\tau_T(A) \ge c_1 \sum_{t=2}^T \|Z_t^* - Z_{t-1}^*\|^2 \ge c_2T$, for some numerical constants $c_1$ and $c_2$. This even holds under strong assumptions on the operators $F_t$, such as strong monotonicity. For a formal result, see Lemma A.1 in Appendix A.

This observation motivates the introduction of a quantity we will use repeatedly to measure the hardness of time-varying VI problems: the quadratic path length of the solutions:

$$ P_T^* = \sum_{t=2}^{T} \|Z_t^* - Z_{t-1}^*\|^2. $$

Since general results for arbitrary sequences are not achievable, we restrict our attention to two important special cases defined below: tame and periodic time-varying VI.

**Definition 1.1** (Tame Time-Varying VI problem). We call a time-varying VI problem tame if the quadratic path length of the solutions $P_T^*$ is sublinear. That is, if there exists time-independent constants $\alpha \in [0, 1)$ and $c > 0$ such that $P_T^* \le cT^\alpha$ for all $T \ge 2$.

For tame time-varying VI problems, we show tracking guarantees for the class of algorithms satisfying a one-step contraction property, in Section 2. In Sections 3 and 4, we focus on periodic time-varying problems, which we define naturally as follows. Note that periodic problems cannot be tame unless the equilibrium is constant.

**Definition 1.2** (Time-varying VI problem with periodic solutions). A time-varying VI problem VIP$((F_t), Z)$ has periodic solutions if there exists $k \in \mathbb{N}$, such that $Z_{t+k}^* = Z_t^*$ for all $t \in \mathbb{N}$.

Periodic time-varying VI problems, that is, time-varying VI problems with $F_{t+k} = F_t$ for all $t \in \mathbb{N}$, are time-varying VI problems with periodic solutions; the converse does not necessarily hold.
and Stampacchia [2000, Chapter 7], or Facchinei and Pang [2003] for more details and examples. Beyond these well-known cases, generalized linear models (GLMs) Juditsky and Nemirovsky [2019] also fit within the framework of monotone VI problems, though they are neither convex optimization problems nor concave games but statistical parameter estimation problems (see Example 4.3 for a time-varying illustration). In Section 4, we discuss several examples in detail and the guarantees our results can provide for these application examples.

We study the problem of online tracking of solutions for time-varying VI. More precisely, we consider a sequence of operators $(F_t)_{t \in \mathbb{N}}$ revealed in an online fashion to a learner. At time step $t$, given the past observations, the learner outputs a candidate solution $Z_t$, and then observes the operator $F_t$. Formally, an online algorithm is a sequence of update rules $A_t: Z \times \mathcal{F}^{t-1} \to Z$ which act online: at time step $t$ the output of the algorithm may depend on $F_1, \dots, F_{t-1}, Z_{t-1}$, but not on $(F_s)_{s \ge t}$. An instance of a time-varying VI problem is characterized by the sequence $(F_t)$ together with the domain $Z$; we denote it by VIP$((F_t)_{t \in [T]}, Z)$.

We focus especially on sequences of operators that all admit a unique solution, denoted by $(Z_t^*)$. A sufficient, but not necessary, condition for uniqueness is that the operator be strongly monotone. For VIP$((F_t)_{t \in [T]}, Z)$, we examine the tracking guarantees of algorithms. Denoting the $t$-th output of algorithm $A$ by $Z_t = A_t(Z_{t-1}, (F_s)_{s \in [t-1]})$, the tracking error of $A$ is

$$ \tau_T(A) = \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2. $$

We exhibit conditions under which this tracking error is sublinear, i.e., there exists an $\alpha \in [0, 1)$ and constant $c \ge 0$, such that $\tau_T(A) \le cT^\alpha$. Note that guaranteeing a sublinear tracking error for time-varying VI problems may be impossible for some sequences, especially when solutions vary too much, see Lemma A.1 in Appendix A. This motivates the introduction of the quadratic path length, which measures the hardness of time-varying VI problems

$$ P_T^* = \sum_{t=2}^{T} \|Z_t^* - Z_{t-1}^*\|^2. $$

Since general results for arbitrary sequences are not achievable, we restrict our attention to two important special cases defined below: tame and periodic time-varying VI.

**Definition 1.1** (Tame Time-Varying VI problem). We call a time-varying VI problem tame if the quadratic path length of the solutions $P_T^*$ is sublinear. That is, if there exists time-independent constants $\alpha \in [0, 1)$ and $c > 0$ such that $P_T^* \le cT^\alpha$ for all $T \ge 2$.

For tame time-varying VI problems, we show tracking guarantees for the class of algorithms satisfying a one-step contraction property, in Section 2. In Sections 3 and 5, we focus on periodic time-varying problems, which we define naturally as follows. Note that periodic problems cannot be tame unless the equilibrium is constant.

**Definition 1.2** (Time-varying VI problem with periodic solutions). A time-varying VI problem VIP$((F_t), Z)$ has periodic solutions if there exists $k \in \mathbb{N}$, such that $Z_{t+k}^* = Z_t^*$ for all $t \in \mathbb{N}$.

Periodic time-varying VI problems, that is, time-varying VI problems with $F_{t+k} = F_t$ for all $t \in \mathbb{N}$, are time-varying VI problems with periodic solutions; the converse does not necessarily hold.

**Summary of contributions.** We provide three sets of contributions to the study of tracking in time-varying VI, also summarized in Table 1:

1. Section 2 focuses on tame time-varying VI problems. We provide path-length bounds for algorithms satisfying a one-step contraction property. This condition is satisfied in a variety of settings, with bounded or unbounded domains, strongly monotone and non-strongly-monotone cases; see Section 2.2 for examples.

2. We then shift our attention to periodic and strongly-monotone sequences of operators, with two main results. We first give an algorithm with a logarithmic tracking bound for bounded domains and bounded operators, in Theorem 3.1. In the case of Lipschitz operators on unbounded domains, we design an algorithm with a constant tracking bound in Theorem 3.2. Both algorithms rely on a meta-algorithm aggregating over a collection of algorithms. Hence, they adapt to the a priori unknown period of the sequences and only require an upper bound on that period.

3. In Section 5, we examine gradient descent with a constant step size on a periodic optimization problem. We show that the composition of the periodic operators yields an autonomous (time-invariant) dynamical system. In this autonomous setting, we observe a range of behaviors, including convergence, periodicity, and Li-Yorke chaos, depending on the step size.
**Examples.** Let us illustrate the practical relevance of tame and periodic VI problems with some examples. See also Facchinei and Pang [2003, Chapter I] for many more examples of applications of variational inequalities, which often have natural time-varying extensions.

**Example 1.1** (Repeated kelly auctions on seasonal markets.) We adapt an example from Duvocelle et al. [2023]. Consider a repeated Kelly auction (Kelly et al. [1998]) where $n$ bidders participate in a repeated auction on a limited good. Let $R_t \ge 0$ denote the minimal price, i.e., the entry barrier set at time $t$. For all $i \in [n]$, given a maximum budget $b_t^{(i)}$, and their value for the good, $v_t^{(i)} \ge 0$, player $i$ bids $x_t^{(i)} \in [0, b_t^{(i)}]$. Each player gets a fraction of the good relative to their bid. Their utility function is

$$ \nu_t^{(i)}(x_t^{(i)}, x_t^{(-i)}) = v_t^{(i)} \frac{x_t^{(i)}}{R_t + \sum_{j \in [n]} x_t^{(j)}} - x_t^{(i)}. $$

For many goods, the demand or the supply depends on exogenous factors that vary with time. For example, the availability of agricultural products depends on the season, and the electricity demand depends on the time of the day. These conditions naturally affect the minimal price $R_t$, and the individual values $v_t^{(i)}$, making the game time-varying with seasonal shifts.

**Example 1.2** (Linear regression on streaming data.) Consider a linear regression problem with a continuous stream of i.i.d. data as introduced in Foster [1991]. That is, let $\{(a_1, b_1), \dots, (a_{n_t}, b_{n_t})\} \subset \mathbb{R}^d \times \mathbb{R}$ denote the data at time $t$, and $A_t := [a_1, \dots, a_{n_t}]^\top \in \mathbb{R}^{n_t \times d}$ and $b_t := [b_1, \dots, b_{n_t}]^\top \in \mathbb{R}^{n_t}$ the corresponding data matrix and vector, respectively. Suppose the data points $(a_i, b_i)$ are coming from a stream of measurements of the same experiment. In an experimental setting, it can be beneficial to have a real-time regression on the data stream (see e.g. the software tool [Mat, 2020]). Hence, define $f_t: \mathbb{R}^d \to \mathbb{R}$ with $f_t(x) = \|A_t x - b_t\|^2 + \lambda \|x\|^2$. Under the assumption that the experimental data follows a fixed sub-gaussian distribution, the solution to the empirical regression problem approaches the solution with respect to the true distribution with high probability, see Shalev-Shwartz and Ben-David [2014, Section 13]. This defines an optimization problem with a slowly shifting solution and thus a *tame time-varying* optimization problem.

The following example illustrates that VI problems are more general than optimization or saddle-point problems.

**Example 1.3** (Generalized linear models.) We present an example from Juditsky and Nemirovsky [2019], noting that this work contains further examples and applications. Consider a generalized linear model corrupted by Gaussian noise; i.e., we obtain observations $(a_t, b_t) \in \mathbb{R}^d \times \mathbb{R}$, $t = 1, \dots, T$, that satisfy $b_t = \phi(\langle Z^*, a_t \rangle) + \xi_t$, where $\phi: \mathbb{R} \to \mathbb{R}$ is a continuous and increasing link function, $\xi_t \sim \mathcal{N}(0, \sigma^2)$, and $Z^* \in \mathbb{R}^d$ is an unknown parameter vector. In this i.i.d. scenario, the maximum likelihood estimator of the unknown parameter $Z \in \mathbb{R}^d$ is given by the (nonlinear) least squares problem

$$ \min_{Z \in \mathbb{R}^d} \left\{ \frac{1}{T} \sum_{t=1}^T [b_t - \phi(\langle Z, a_t \rangle)]^2 \right\}. $$

Note this problem is typically nonconvex (with the exception in the case the link function is affine). Juditsky and Nemirovsky [2019] observe that despite the nonconvexity above, the parameter estimation problem can be reformulated as a monotone VI problem VIP($F, \mathbb{R}^d$) with operator

$$ F(Z) = \frac{1}{T} \sum_{t=1}^{T} a_t [\phi(\langle Z, a_t \rangle) - b_t]. \quad (1) $$

Moreover, under minimal distributional assumptions it can be proved that the population-type operator (i.e., the expected value of $F$ under the target distribution) becomes strongly monotone. The above VI can also be considered in a streaming setting, such as proposed in Example 1.2. If needed, regularization can be introduced by adding $\lambda Z$ to the operator in equation (1).

**On naive solutions to time-varying problems.** A first approach one should consider is to simply ignore the time-varying nature of the problem. Indeed for fluctuations small enough, the effect of the time variations on the location of the solutions could be negligible. More critically, detecting the variations can be a difficult task in itself. In games, the variations of the loss of a player could be due to fluctuations in their losses, or to a change of other player's actions. Noisy feedback could be another source of uncertainty regarding the time variations of the objective.

For these reasons, we examine some examples of algorithms oblivious to the periodicity of the sequence. Let us illustrate one elementary yet striking consequence of ignoring the time fluctuations, on a sequence of 1-dimensional
<table><caption>Table 1: Overview of Results. Everywhere, $c, c'$ denote constants independent of the number of iterations $T$. The constant $D_0 = \max_{i \in [k]} \|Z_1 - Z_i^*\|$ denotes the maximal distance of the initial iterate to the solutions $Z_1^*, \dots, Z_k^*$, and $K$ denotes a known upper bound on $k$. $D_0$ and $K$ are instance dependent constants independent of $T$. The abbreviations of the algorithms are: EW: Exponential Weights, OGD: Online Gradient descent, GD: Gradient Descent. For details on the algorithms, see Section 3.</caption><thead><tr><th></th><th>Algorithm A</th><th>Results</th></tr></thead><tbody><tr><td rowspan="2">Tame VIP</td><td>Any contractive algorithm, e.g.:<br/>- Gradient descent<br/>- Resolvent iteration</td><td><b>Theorem 2.1:</b> Sublinear tracking error<br/>$\tau_T(A) \le cT^\alpha$.<br/><b>Assumptions:</b> tame (Def. 1.1), contractive (Def. 2.1). <b>Remark:</b> strong monotonicity not necessarily required (cf. Example 2.3).</td></tr><tr><td>Meta-algorithm with:<br/>Aggregation: EW, constant learning rate<br/>Base Algorithm: OGD</td><td><b>Theorem 3.1:</b> Logarithmic tracking error<br/>$\tau_T(A) \le c \frac{k(G+D)^2}{\mu^2} \log T + c'$.<br/><b>Assumptions:</b> k-periodic (Def. 1.2), bounded domain, strongly monotone ($\mu$-Mono) and bounded operators<br/><b>Remark:</b> Computational cost per iteration is constant.</td></tr><tr><td rowspan="3">k-Periodic VIP</td><td>Meta-algorithm with:<br/>Aggregation: EW, adaptive learning rate<br/>Base Algorithm: GD</td><td><b>Theorem 3.2:</b> Constant tracking error<br/>$\tau_T(A) \le c k D_0^2 \frac{L^4}{\mu^4} \log K + c'$.<br/><b>Assumptions:</b> k-periodic (Def. 1.2), L-Lipschitz, $\mu$-strongly monotone operators.<br/><b>Remark:</b> Domain can be unbounded. Computational cost per iteration is in the order of $K$.</td></tr><tr><td>GD, constant step-size</td><td>Phase I: Convergence.<br/>Phase II: Composite phase with convergence, periodicity, chaos, and divergence. (Obs. 1)<br/>Phase III: Divergence.</td></tr></tbody></table>

We stress that the emergence of chaos, as discussed in Section 5, is a property of the autonomous dynamical system induced by the periodic composition of operators, and not of general non-autonomous (time-varying) problems. We illustrate the practical consequences of our tracking guarantees with examples in Section 4.

**Remark 1.1** (Ignoring Time-Variations and Practical Relevance). In the context of time-varying problems with uncertainty about their structure, a first (arguably reasonable) attempt would be to ignore the problem's time-varying nature. If fluctuations are small, their impact on solutions may be negligible.

With this in mind, we examine a basic example of gradient descent with a fixed step-size applied on a periodic sequence of strongly convex functions with a shared minimizer (see Appendix D), and show that tuning the learning rate for convergence is hard despite the simplicity of the example. Section 5 furthers this investigation, showing chaotic phenomena even in the case where all periodic functions share the same minimizer. These negative results in such simple settings justify and motivate the more involved approach carried out in Section 3.

## 1.1 Dynamic regret minimization.

Note that showing tracking bounds for an algorithm is related to *dynamic regret minimization* [Zinkevich, 2003] in the special case of time-varying optimization problems. Indeed, consider a sequence of functions $f_1, \dots, f_T$ and a sequence of comparators $u_1, \dots, u_T \in \mathbb{Z}$. These comparators can be the minimizers $u_t^*$ of the functions $f_t$. The dynamic regret of an online algorithm $\mathcal{A}$ is

$$d\mathrm{Reg}_T^A((u_t)_{t \in [T]}) := \sum_{t=1}^T (f_t(x_t) - f_t(u_t)).$$
<table><caption>Table 1: Overview of Results. Everywhere c, c' denote constants independent of the number of iterations T. The constant $D_0 = \max_{i \in [k]} \|Z_1 - Z_i^*\|$ denotes the maximal distance of the initial iterate to the solutions $Z_1^*, \dots, Z_k^*$, and K denotes a known upper bound on k. $D_0$ and K are instant dependent constants independent of T. The abbreviations of the algorithms are: EW: exponential weights, OGD: online gradient descent, GD: gradient descent. For details on the algorithms, see Section 3.</caption><thead><tr><th></th><th>Algorithm A</th><th>Results</th></tr></thead><tbody><tr><td rowspan="2">Tame Time-Varying VIP</td><td>Contractive algorithms. E.g.<br/>- Gradient descent<br/>- Resolvent iteration</td><td>Theorem 2.1: Sublinear tracking error<br/>τ<sub>T</sub>(A) ≤ cT<sup>α</sup>.<br/>Assumptions: tame (Def. 1.1), contractive (Def. 2.1)</td></tr><tr><td>Meta-algorithm with:<br/>Aggregation: EW, constant<br/>learning rate<br/>Base Algorithm: OGD</td><td>Theorem 3.1: Logarithmic tracking error<br/>τ<sub>T</sub>(A) ≤ ck<sup>1</sup>/μ<sup>2</sup> log T + c'.<br/>Assumptions: k-periodic (Def. 1.2), bounded domain, strongly<br/>monotone (μ-Mono)</td></tr><tr><td rowspan="2">k-Periodic VIP</td><td>Meta-algorithm with:<br/>Aggregation: EW, adaptive<br/>learning rate<br/>Base Algorithm: GD</td><td>Theorem 3.2: Constant tracking error<br/>τ<sub>T</sub>(A) ≤ ckD<sub>0</sub><sup>2</sup>L<sup>4</sup>/μ<sup>4</sup> log K + c'.<br/>Assumptions: k-periodic (Def. 1.2), unconstrained, L-Lipschitz<br/>(L-Lip), strongly monotone (μ-Mono)</td></tr><tr><td>GD with constant step-size</td><td>Phase I : Convergence.<br/>Phase II : Composite phase with convergence, periodic,<br/>chaotic, and divergence. (Obs. 1)<br/>Phase III: Divergence.</td></tr></tbody></table>

quadratic problems, using the gradient descent algorithm with a fixed step size. We will see that increasing the learning rate can move the algorithm from a diverging regime to a converging one, and even to the fastest convergence rates. This goes against the intuitive idea that decreasing the learning should make the sequence more stable.

For $\mathcal{Z} = \mathbb{R}$, consider the sequence of functions $(f_t)$ defined by $f_{2t}(x) = x^2/2$ and $f_{2t+1}(x) = 4x^2$, which corresponds to the periodic time-varying VIP((F_t), $\mathbb{R}$) with $F_{2t}(x) = x$ and $F_{2t+1}(x) = 8x$, and with a constant solution $x^* \equiv 0$. We observe that $(F_t)$ defines a 2-periodic VI problem and a VI problem with a fixed (i.e., 1-periodic) solution. Both of the operators $F_{2t}$, $F_{2t+1}$ are 1-strongly monotone, and 8-Lipschitz. Consider gradient descent with step size $\eta$, which consists in selecting at all times $t$ the sequence defined by $x_{t+1} = x_t - \eta F_t(x_t)$. The actions at every odd time satisfy

$$x_{t+2} = (1 - \eta)(1 - 8\eta)x_t.$$

Gradient descent converges to the optimum 0 (for generic initial points) if and only if the common ratio of the geometric sequence is less than 1, i.e., if $|1 - \eta| |1 - 8\eta| < 1$. Moreover, when the sequence converges, the common ratio dictates its speed of convergence to the optimum: the closer it is to 0, the faster the convergence. In particular, when tuned with $\eta = 1$, gradient descent converges to the optimum in a single step. On the other hand, there is a range of smaller step sizes for which the sequence diverges. For example, for $\eta = 1/2$, the common ratio is $3/2 > 1$, and the iterates diverge.

We extend this example and describe much richer phenomena on examples of operators beyond the linear case, showing that the behavior of an algorithm with static guarantees may be wildly unpredictable on simple periodic time-varying problems.
Note that the comparators ($u_t$) are not necessarily the minimizers of the functions ($f_t$). Guarantees for dynamic regret where the comparators are restricted to the minimizers $u_t^*$ of strongly convex and smooth functions ($f_t$) have been studied, e.g. in Zhao and Zhang [2021]. Relating these results to tracking, we note that the quantity $(2/\mu)$ dReg$_T^A((u_t^*)_t \in [T])$ is an upper bound on the tracking error $\tau_T(\mathcal{A})$ over the set of $\mu$-strongly convex functions $\hat{F}_\mu$. Hence, if an algorithm has sublinear dynamic regret, it is a tracking algorithm. Dynamic regret bounds typically depend on the path length $P_T^1 = \sum_{t=1}^{T-1} \|u_t - u_{t+1}\|$, that is, the norm differences of the comparator terms. Indeed, existing lower bounds on the dynamic regret show that $\sqrt{T(P_T^1 + 1)}$ is unavoidable for convex functions Zhang et al. [2018]. Hence, showing sublinear bounds require additional assumptions, either on the path length $P_T^1$ or on the distances between the functions $f_t$.

## 1.2 Related work.

Tracking equilibria in time-varying games has received growing interest in the game theory community in recent years, see Rivera Cardoso et al. [2019], Fiez et al. [2021], Mertikopoulos and Staudigl [2021], Duvocelle et al. [2023], or Feng et al. [2023]. Duvocelle et al. [2023, Theorem 2 and Corollary 3] provided guarantees for tracking equilibria in time-varying games under a very general feedback model that allows for noisy and biased observations. In the tame case, they provide tracking bounds for a proximal point method under strong monotonicity. Yan et al. [2023] strengthen these results in the exact feedback model under an additional smoothness assumption, providing tighter and more adaptive bounds. Feng et al. [2023] point out a surprising example in time-varying bilinear games, in which the extra-gradient method outperforms optimistic methods of Rakhlin and Sridharan [2013]; these two algorithms were often thought to enjoy almost identical behavior.

A related line of work studies no-regret dynamics in games, focusing on algorithms designed to minimize regret at the individual player level. In repeated static strongly-convex-strongly-concave zero-sum games, the sum of the regret of the players bounds the distance of the joint average plays to the equilibrium of the game (see, e.g., Example 3.2). Zhang et al. [2022] investigates the impact of time variations of the game, and derives tracking guarantees in the tame case when players use dynamic regret algorithms.

As seen in the previous section, tracking solutions in the special case of optimization are related to dynamic regret minimization. The literature distinguishes dynamic regret against arbitrary comparators (Zinkevich [2003], Baby and Wang [2022]) and against the minimizers of the sequences of functions observed (Besbes et al. [2015]); the former is a more demanding objective, while the latter suffices for tracking when functions are strongly convex. In this case, the dynamic regret against minimizers bounds the tracking error up to a multiplicative factor. For strongly convex and smooth functions, Mokhtari et al. [2016], using the fact that gradient descent is contractive (see Definition 2.1), provides tracking bounds of order $P_T^*$. Our analysis in Theorem 2.1 is a generalization of their results to VI problems and to other algorithms. For an in-depth comparison, see Section 2.4.

Online periodic optimizations and periodic games are less explored topics. Fiez et al. [2021] show that the continuous-time dynamics of gradient-descent-ascent are Poincaré recurrent in periodic bilinear unconstrained zero-sum games, showing a form of stability despite the time variations—note that the equilibrium is constant and equal to 0 in their analysis. In that same setting, but with discrete time, Feng et al. [2023] show the convergence of extra-gradient methods to the constant equilibrium. Our results are complementary: we consider the strongly monotone setting, and provide tracking bounds.

The dynamics of optimization methods with large step sizes have received attention recently, with some works exhibiting chaotic behavior. Chaos may also appear in natural game dynamics, e.g. in Piliouras and Yu [2023], Falniowski and Mertikopoulos [2024]. Closest to our results of Section 5 is that of Chen et al. [2023], which characterize different phases in the behavior of gradient descent in quadratic regression, as the learning rate increases. To the best of our knowledge, we are the first to exhibit such a phenomenon in the periodic time-varying case. We also note that periodic time-varying VI problems are related to incremental gradient descent Nedic and Bertsekas [2001]. Our observations might also be of interest in this context.

Furthermore, two related lines of research in the area of game theory are worth noting. Mechanism design [Baudin et al., 2023, Deng et al., 2021] and often specifically auction design [Nedelec et al., 2022, Rahme et al., 2021] frameworks, in which a central entity modifies the parameters of a game, induce time-varying games for the players. Stochastic games model repeated settings in which the joint actions of players influence future payoffs, causing time variations in the individual objectives. Our results may formally apply to some specific settings in these models, but they ignore the crucial strategic aspects involved in planning for the long-term consequences of current actions. In Section 5.3, we also mention a relation to the theory of Iterated Function Systems (IFS), see e.g., Falconer [1990, Chapter 9]. Relations between stochastic gradient descent methods and IFS were recently used in Hodgkinson et al. [2021] to es-
## Summary of contributions.

1. Section 2 focuses on tame time-varying VI problems. In Theorem 2.1, we establish a sufficient condition for tracking, namely that the algorithm used is *contractive* (see Definition 2.1). This generalizes an argument by Mokhtari et al. [2016] from online optimization to VI problems. In Section 2.2 we provide examples of contractive algorithms beyond the strongly monotone and Lipschitz case, showing that the criterion of having a contractive algorithm reaches beyond several existing settings.

2. For periodic VI problems, results from Section 2 are vacuous since the path length is always linear. In Section 3, Theorem 3.1 we give logarithmic tracking bounds for strongly-monotone periodic time-varying VI problems on a bounded domain. Under the additional assumption that the operators are Lipschitz, Theorem 3.2 improves these results for periodic time-varying VI problems to constant bounds for the unconstrained case. Our new algorithm carefully combines online learning techniques—surrogates losses and aggregation [Vovk, 1995, Van Erven et al., 2021]—to adapt to the unknown period of the sequence via a reduction from VI problems to online convex optimization.

3. Contrasting these results, in Section 4, we show that the dynamical system induced by gradient descent with fixed step size and a time-varying VI problems of period two can exhibit chaotic behavior. This shows the necessity of designing specific algorithms for the periodic setting. In an experimental study, we demonstrate that depending on the step size, several phases may occur: (1) With a sufficiently small step size we observe convergence. (2) For a certain interval of step-sizes, we observe chaotic, periodic, divergence, or convergence behavior. (3) For a sufficiently large step size, the system diverges. The occurrence of chaos per se may not be surprising to some readers: similar observations were made for example in the context of gradient descent Chen et al. [2023]. However, in existing work, a stringent transition from convergence to periodic to chaos to divergence was observed. We observe that for time-varying problems, this phase landscape is richer. For example, increasing the step-size leads to convergence again after observing chaotic behavior for a smaller step size.

## 1.1 Dynamic regret minimization and equilibrium tracking in time-varying games.

Note that showing tracking bounds for an algorithm is related to *dynamic regret minimization* [Zinkevich, 2003] in the special case of time-varying optimization problems. Indeed, consider a sequence of functions $f_1, \dots, f_T$ and a sequence of comparators $u_1, \dots, u_T \in \mathcal{Z}$. These comparators can be the minimizers $u_t^*$ of the functions $f_t$. The dynamic regret of an online algorithm $\mathcal{A}$ is

$$d\mathrm{Reg}_T^A ((u_t)_{t \in [T]}) := \sum_{t=1}^T f_t(x_t) - f_t(u_t).$$

The quantity $(2/\mu)$ $d\mathrm{Reg}_T^A (\{u_t^*\}_{t \in [T]})$ is an upper bound on the tracking error $\tau_T(\mathcal{A})$ over the set of $\mu$-strongly convex functions $\hat{\mathcal{F}}_\mu$. Hence, if an algorithm has sublinear dynamic regret, it is a tracking algorithm. Dynamic regret bounds typically depend on the path length $P_T^1 = \sum_{t=1}^{T-1} \|u_t - u_{t+1}\|$, that is, the norm differences of the comparator terms. Indeed, existing lower bounds on the dynamic regret show that $\sqrt{T(P_T^1 + 1)}$ is unavoidable for convex functions Zhang et al. [2018]. Hence, showing sublinear bounds requires additional assumptions, either on the path length $P_T^1$ or on the distances between the functions $f_t$.

Note, that it is not necessarily required that the comparators $(u_t)$ are the minimizers of the functions $(f_t)$. Indeed, the guarantees for dynamic regret where the comparators are restricted to the minimizers $u_t^*$ of strongly convex and smooth functions $(f_t)$ manifest stronger problem dependencies Zhao and Zhang [2021].

Next, consider the special case of time-varying games and the tracking of equilibria. This problem received a growing interest in the game theory community in recent years, see Rivera Cardoso et al. [2019], Fiez et al. [2021], Duvocelle et al. [2023], or Feng et al. [2023]. To establish tracking guarantees for the Nash-equilibria of time-varying games, Duvocelle et al. [2023] assumed that the games are *strongly monotone*: That is, that the associated operator given by the pseudo-gradient of the utility functions is strongly monotone. Under the assumption that a game satisfies $(\mu$-Mono), bounds on the tracking error $\tau_T(\mathcal{A})$ are shown in Duvocelle et al. [2023]. Under the additional assumption that the variations of the game are sublinear, these bounds give tracking guarantees for various algorithms. We note that $\mu$-Mono is implied if the utility functions $\nu^{(i)}(\cdot, x^{(-i)})$ are strongly convex for all players $i$ are with respect to any $x^{(-i)} \in \mathcal{Z}^{(-i)}$. In the special case of a two-player zero-sum game (a.k.a. a saddle-point problem), we refer to this condition as *strongly-convex-strongly-concave*.
tablish convergence or generalization guarantees for popular machine learning methods. We build on similar relations to establish a link between existing results for IFS and our experimental observations.

## 2 Tracking by contractive algorithms.

We derive tracking bounds for time-varying variational inequalities using a contraction property of algorithms. While this property holds for many common first-order algorithms for strongly monotone VIs or VIs with bounded sets $Z$, we also show that strong monotonicity or boundedness is not necessarily required for contractiveness. In Section 2.1, we define contractive algorithms and show that they automatically enjoy tracking bounds (Theorem 2.1). Examples of such algorithms are presented in Section 2.2, followed by a discussion of their implications in Section 2.4.

### 2.1 Contractive implies tracking.

In strongly convex and smooth optimization, the gradient descent updates ensure that the distance to the minimizer decreases by a constant multiplicative factor at every time step. Mokhtari et al. [2016] use this property to derive dynamic regret bounds for gradient descent in online optimization. We generalize this observation to time-varying VI problem, and contractive algorithm, which are defined as follows:

**Definition 2.1 (Contractive Algorithms).** Let $C \in (0, 1)$ and $A$ be an algorithm composed of update rules ($A_t$). The update rule $A_t$ is said to be C-contractive over a set of operators $\mathcal{F}$ if for all $F \in \mathcal{F}$ with solution $Z^* \in Z$, for any $Z \in Z$,

$$ \|Z_+ - Z^*\| \le C \|Z - Z^*\|, $$

where we denoted by $Z_+$ the image of $Z$ under the update rule $A_t$ after observing $F$. We say that the algorithm $A$ is C-contractive over $\mathcal{F}$ if $A_t$ is C-contractive for all $t \ge 1$.

**Theorem 2.1.** Suppose $A$ is C-contractive over $\mathcal{F}$. Then for any sequence of operators in $\mathcal{F}$ with solutions $(Z_t^*)_{t \in [T]}$, the tracking error is bounded by

$$ \tau_T(\mathcal{A}) \le \frac{1}{(1-C)^2} \sum_{t=2}^{T} \|Z_t^* - Z_{t-1}^*\|^2 + \frac{1}{1-C} \|Z_1 - Z_1^*\|^2. $$

*Proof.* Proof. Using the inequality $\|u + v\|^2 \le (1 + \alpha)\|u\|^2 + (1 + \alpha^{-1})\|v\|^2$ for any $u, v \in \mathbb{R}^d$ and $\alpha > 0$,

$$ \begin{aligned} \|Z_{t+1} - Z_{t+1}^*\|^2 &\le \frac{1}{C} \|Z_{t+1} - Z_t^*\|^2 + \left(1 + \frac{1}{1/C - 1}\right) \|Z_{t+1}^* - Z_t^*\|^2 \\ &\le C \|Z_t - Z_t^*\|^2 + \frac{1}{1-C} \|Z_{t+1}^* - Z_t^*\|^2, \end{aligned} $$

where the second inequality follows by the contraction property. Summing over $t \in [T-1]$ yields

$$ \sum_{t=1}^{T-1} \|Z_{t+1} - Z_{t+1}^*\|^2 \le C \sum_{t=1}^{T-1} \|Z_t - Z_t^*\|^2 + \frac{1}{1-C} \sum_{t=1}^{T-1} \|Z_{t+1}^* - Z_t^*\|^2. $$

Then, adding $\|Z_1 - Z_1^*\|^2$ to both sides, and upper bounding the first sum on the right-hand side,

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le C \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 + \frac{1}{1-C} \sum_{t=1}^{T-1} \|Z_{t+1}^* - Z_t^*\|^2 + \|Z_1 - Z_1^*\|^2. $$

Therefore, after reorganizing the sum,

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le \frac{1}{(1-C)^2} \sum_{t=1}^{T-1} \|Z_{t+1}^* - Z_t^*\|^2 + \frac{1}{1-C} \|Z_1 - Z_1^*\|^2, $$

which completes the proof.

<figure><img src="image_1.png" /></figure>

<figure><img src="image_2.png" /></figure>

Note that Theorem 2.1 does not require $Z$ to be bounded. However, some algorithms require a bounded convex set $Z$ in order to satisfy the contraction property (Definition 2.1), e.g., Example 2.4.
## 1.2 Related work.

Duvocelle et al. [2023, Theorem 2 and Corollary 3] introduce the problem of tracking equilibria in time-varying games under a very general feedback model that allows for noisy and biased observations. In the same case, they provide tracking bounds for a proximal point method under strong monotonicity. Yan et al. [2023] strengthen these results in the exact feedback model under an additional smoothness assumption, providing tighter and more adaptive bounds. Feng et al. [2023] point out a surprising example in time-varying bilinear games, in which the extra-gradient method outperforms optimistic methods of Rakhlin and Sridharan [2013]; these two algorithms were often thought to enjoy almost identical behavior.

A related line of work studies no-regret dynamics in games, focusing on algorithms designed to minimize regret at the individual player level. In repeated static strongly-convex-strongly-concave zero-sum games, the sum of the regret of the players bounds the distance of the joint average plays to the equilibrium of the game (see, e.g., Example 3.2). Zhang et al. [2022] and Anagnostides et al. [2024] investigate the impact of time variations of the game, and derive tracking guarantees in the same case when players use dynamic regret algorithm. See Section 2.4 for a detailed discussion on the rates achieved in these different models and a comparison with ours.

Tracking solutions in the special case of optimization is also related to dynamic regret minimization. Literature distinguishes dynamic regret against arbitrary comparators (Zinkevich [2003], Baby and Wang [2022]) and against the minimizers of the sequences of functions observed (Besbes et al. [2015]); the former is a more demanding objective, but the latter sometimes suffices for tracking. Indeed, for strongly convex functions, which always have a unique minimizer, the dynamic regret against the sequence of minimizers upper bounds the tracking error up to a multiplicative factor. For strongly convex and smooth functions, Mokhtari et al. [2016], using the fact that gradient descent is contractive (see Definition 2.1), provides tracking bounds of order $P_T^*$. Our analysis in Theorem 2.1 is a generalization of their results to VI problems and to other algorithms.

Online periodic optimizations and periodic games are less explored topics. Fiez et al. [2021] show that the continuous-time dynamics of gradient-descent-ascent are Poincaré recurrent in periodic bilinear unconstrained zero-sum games, showing a form of stability despite the time variations—note that the equilibrium is constant and equal to 0 in their analysis. In that same setting, but with discrete time, Feng et al. [2023] show the convergence of extra-gradient methods to the constant equilibrium. Our results are complementary: we consider the strongly monotone setting, and provide tracking bounds.

The dynamics of optimization methods with large step sizes have received attention recently, with some works exhibiting chaotic behavior. Chaos may also appear in natural game dynamics, e.g. in Piliouras and Yu [2023], Falniowski and Mertikopoulos [2024]. Closest to our results of Section 4 is that of Chen et al. [2023], which characterize different phases in the behavior of gradient descent in quadratic regression, as the learning rate increases. To the best of our knowledge, we are the first to exhibit such a phenomenon in the periodic time-varying case. We also note that periodic time-varying VI problems are related to single-shuffle SGD Nguyen et al. [2021], Mishchenko et al. [2020] as well as to incremental gradient descent Nedic and Bertsekas [2001]. Hence, our observations might also be of interest in this context.

In Section 4.3, we also mention a relation to the theory of Iterated Function Systems (IFS), see e.g., Falconer [1990, Chapter 9]. Relations between stochastic gradient descent methods and IFS were recently used in Hodgkinson et al. [2021] to establish convergence or generalization guarantees for popular machine learning methods. We build on similar relations to establish a link between existing results for IFS and our experimental observations.

# 2 Tracking by contractive algorithms.

We show how to derive tracking bounds in time-varying variational inequalities, from a contraction property satisfied by common algorithms on static problems with assumptions like strong monotonicity. In Section 2.1, we define *contractive* algorithms and show that they automatically enjoy tracking bounds (Theorem 2.1). We then provide examples of contractive algorithms in Section 2.2, and discuss consequences in Section 2.4.

## 2.1 Contractive implies tracking.

In strongly convex and smooth optimization, the gradient descent updates ensure that the distance to the minimizer decreases by a constant multiplicative factor at every time step. Mokhtari et al. [2016] use this property to derive dynamic regret bounds for gradient descent in online optimization. We generalize this observation to time-varying VI problem, and contractive algorithm, which are defined as follows:
## 2.2 Examples of contractive algorithms.

The contraction property (Definition 2.1) is satisfied in a wide range of settings, often in the presence of curvature of the objective (in the form of strong monotonicity of the operator), but not always. Specifically, Example 2.3 shows that strong monotonicity is not a necessary assumption for contractivity. Furthermore, note that some examples do not require a bounded $\mathcal{Z}$.

**Example 2.1.** Let $F: \mathcal{Z} \to \mathbb{R}^d$ be a $\mu$-strongly-monotone and $L$-Lipschitz operator, where $\mathcal{Z} \subseteq \mathbb{R}^d$ is a closed convex set. Let $T_\eta: z \mapsto \text{Proj}_\mathcal{Z}(z - \eta F(z))$ denote the projected forward operator for some step-size $\eta > 0$. Then $Z^*$ is a solution to VIP($F$, $\mathcal{Z}$) if and only if it is a fixed point of $T_\eta$ (see Bauschke and Combettes [2011, Proposition 23.38]). Moreover, for $\eta = \mu/L^2$, the operator $T_\eta$ is a contraction: for any $Z_1, Z_2 \in \mathcal{Z}$, using the nonexpansiveness of Euclidean projections, strong-monotonicity and the Lipschitz property of $F$,

$$ \begin{aligned} \|T_\eta(Z_1) - T_\eta(Z_2)\|^2 &\le \|Z_1 - \eta F(Z_1) - (Z_2 - \eta F(Z_2))\|^2 \\ &= \|Z_1 - Z_2\|^2 - 2\eta \langle Z_1 - Z_2, F(Z_1) - F(Z_2) \rangle + \eta^2 \|F(Z_1) - F(Z_2)\|^2 \\ &\le (1 + \eta^2 L^2 - 2\eta\mu) \|Z_1 - Z_2\|^2 = (1 - \kappa^{-2}) \|Z_1 - Z_2\|^2, \end{aligned} $$

where $\kappa = \mu/L$ denotes the condition number. This implies that the projected forward algorithm with step-size $\eta = \mu/L^2$ is a contractive algorithm. In particular, this includes gradient descent for optimization and gradient-descent-ascent for saddle-point optimization.

**Example 2.2.** Suppose $F: \mathbb{R}^d \to \mathbb{R}^d$ is a maximally $\mu$-strongly monotone operator. Let $R$ denote the resolvent operator, that is, $R = (\text{Id} + F)^{-1}$. (Note that $R$ is related to the backward iteration, and that operator $R$ coincides with the proximal mapping when $F$ is a gradient operator.) Then a point $Z^*$ is a solution to the VIP($F$, $\mathbb{R}^d$) if and only if $F(Z^*) = 0$, or equivalently if and only if $R(Z^*) = Z^*$. By Theorem 2.1. (xi) in Bauschke et al. [2011] the resolvent operator is $1/(1+\mu)$-contractive This implies in particular that for any $Z \in \mathbb{R}^d$,

$$ \|R(Z) - Z^*\| = \|R(Z) - R(Z^*)\| \le (1 + \mu)^{-1} \|Z - Z^*\|. $$

**Example 2.3.** Consider the following zero-sum game:

$$ \min_{y \in \mathbb{R}} \max_{x \in \mathbb{R}} \left\{ x^2 + 3 \sin^2(x) + a \sin^2(x) \sin^2(y) - y^2 - 3 \sin^2(y) \right\}, \quad (1) $$

where we assume $a \in [0, 1]$. Note that the game in (1) is not a strongly monotone game and it does not satisfy the diagonal concavity assumption of Rosen [1965] (see Proposition B.2 in Appendix B.1 for details). However, note that for any $y \in \mathbb{R}$, the function $\nu(\cdot, y) = x^2 + 3 \sin^2(x) + a \sin^2(x) \sin^2(y) - y^2 - 3 \sin^2(y)$ satisfies the so called restricted secant inequality (RSI) (see Proposition B.1, Appendix B.1 for details and definition). We show in Proposition B.3 in Appendix B.1, that contraction is satisfied for gradient descent-ascent. A sequence $(a_t)$, with $a_t \in [0, 1]$ defines a time-varying instance of this zero-sum game. Note that this example is inspired and heavily builds on well-known results in optimization by Karimi et al. [2016].

**Example 2.4.** We conclude these examples with a more informal one, regarding the extra-gradient method on bi-affine saddle-point problems, which include in particular equilibria of matrix zero-sum games. Suppose that $\mathcal{Z}$ is the product of two polytopes $\mathcal{P}_1, \mathcal{P}_2$ and that $F$ is a map of the form $F(x, y) = c^\top x + d^\top y - x^\top A y$; that is, $F$ is bi-affine. Korpelevich [1976] introduced the extra-gradient method for saddle-point optimization, which consists in playing

$$ Z_{t+1} = \text{Proj}_{\mathcal{Z}} \{ Z_t - \eta F(\text{Proj}_{\mathcal{Z}} \{ Z_t - \eta F(Z_t) \}) \}. $$

She showed (in Theorem 3 in the reference) that if $F$ has a unique solution $Z_F^*$, then the extra-gradient method is eventually a contraction. Precisely, there exists a (relative) neighborhood of the solution, denoted by $\mathcal{U}(Z_F^*)$, such that extra-gradient becomes a contraction as soon as it reaches $\mathcal{U}(Z_F^*)$. Moreover, suppose in addition that the matrix $A$ is a square matrix and the equilibrium lies in the relative interior of both $\mathcal{P}_1$ and $\mathcal{P}_2$, i.e., in relint($\mathcal{P}_1$) $\times$ relint($\mathcal{P}_2$). In that case, the contraction constant depends only on upper bounds of the euclidean operator norms of $A$, and of $A^{-1}$. To fully characterize the contraction property and complete this example, one would need to explicitly describe $\mathcal{U}(Z_F^*)$. Then, given a subset $\mathcal{C} \subset \mathcal{Z}$, the extra-gradient method started at $Z_0 \in \mathcal{C}$ is a contractive algorithm on the set of VI problems VIP($F$, $\mathcal{Z}$) such that $\mathcal{U}(Z_F^*) \subset \mathcal{C}$. Giving an exact description of $\mathcal{U}(Z_F^*)$ is technically challenging and goes beyond the scope of this article. Recent works (e.g., Wei et al. [2021]) provide finite-time linear convergence bounds with explicit constants, but use arguments that do not explicitly rely upon the contraction property.
**Definition 2.1** (Contractive Algorithms). Let $C \in (0, 1)$ and $A$ be an algorithm composed of update rules ($A_t$). The update rule $A_t$ is said to be C-contractive over a set of operators $\mathcal{F}$ if for all $F \in \mathcal{F}$ with solution $Z^* \in \mathcal{Z}$, for any $Z \in \mathcal{Z}$,

$$ \|Z_+ - Z^*\| \le C \|Z - Z^*\|, $$

where we denoted by $Z_+$ the image of $Z$ under the update rule $A_t$ after observing $F$. We say that the algorithm $A$ is C-contractive over $\mathcal{F}$ if $A_t$ is C-contractive for all $t \ge 1$.

**Theorem 2.1.** Suppose $A$ is C-contractive over $\mathcal{F}$. Then for any sequence of operators in $\mathcal{F}$ with solutions $(Z_t^*)_{t \in [T]}$, the tracking error is bounded by

$$ \tau_T(A) \le \frac{1}{(1-C)^2} \sum_{t=2}^{T} \|Z_t^* - Z_{t-1}^*\|^2 + \frac{1}{1-C} \|Z_1 - Z_1^*\|^2. $$

*Proof.* Using the inequality $\|u+v\|^2 \le (1+\alpha)\|u\|^2 + (1+\alpha^{-1})\|v\|^2$ for any $u, v \in \mathbb{R}^d$ and $\alpha > 0$,

$$ \begin{aligned} \|Z_{t+1} - Z_{t+1}^*\|^2 &\le \frac{1}{C} \|Z_{t+1} - Z_t^*\|^2 + \left(1 + \frac{1}{1/C - 1}\right) \|Z_{t+1}^* - Z_t^*\|^2 \\ &\le C \|Z_t - Z_t^*\|^2 + \frac{1}{1-C} \|Z_{t+1}^* - Z_t^*\|^2, \end{aligned} $$

where the second inequality follows by the contraction property. Summing over $t \in [T-1]$ yields

$$ \sum_{t=1}^{T-1} \|Z_{t+1} - Z_{t+1}^*\|^2 \le C \sum_{t=1}^{T-1} \|Z_t - Z_t^*\|^2 + \frac{1}{1-C} \sum_{t=1}^{T-1} \|Z_{t+1}^* - Z_t^*\|^2. $$

Then, adding $\|Z_1 - Z_1^*\|^2$ to both sides, and upper bounding the first sum on the right-hand side,

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le C \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 + \frac{1}{1-C} \sum_{t=1}^{T-1} \|Z_{t+1}^* - Z_t^*\|^2 + \|Z_1 - Z_1^*\|^2. $$

Therefore, after reorganizing the sum,

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le \frac{1}{(1-C)^2} \sum_{t=1}^{T-1} \|Z_{t+1}^* - Z_t^*\|^2 + \frac{1}{1-C} \|Z_1 - Z_1^*\|^2, $$

which completes the proof.

<figure><img src="image_2.png" /></figure>

**Remark 2.1** (Generality beyond Euclidean vector spaces.) While in the introduction, we defined the VI problem with respect to Euclidean vector spaces for the sake of simplicity, we note that the results from Section 2.1 can be generalized. For example: Let $\mathcal{H}$ denote a real Hilbert space with inner product $\langle \cdot, \cdot \rangle: \mathcal{H} \times \mathcal{H} \to \mathbb{R}$ and induced norm $\| \cdot \|$. Consider the corresponding time-varying VI problem in this generalized space. Theorem 2.1 applies and together with Example 2.2 in the subsequent section, we can obtain tracking guarantees of resolvent operators for tame time-varying strongly monotone VI problem in Hilbert space.

## 2.2 Examples of contractive algorithm.

The contraction property (Definition 2.1) is satisfied in a wide range of settings, often in the presence of curvature of the objective (in the form of strong monotonicity of the operator), but not always. Here are a few standard examples:

**Example 2.1.** Let $F: \mathcal{Z} \to \mathbb{R}^d$ be a $\mu$-strongly-monotone and $L$-Lipschitz operator, where $\mathcal{Z} \subseteq \mathbb{R}^d$ is a closed convex set. Let $T_\eta: z \mapsto \text{Proj}_z(z - \eta F(z))$ denote the projected forward operator for some step-size $\eta > 0$. Then $Z^*$ is a solution to VIP($F$, $\mathcal{Z}$) if and only if it is a fixed point of $T_\eta$ (see Bauschke and Combettes [2011, Proposition 23.38]). Moreover, for $\eta = 1/L$, the operator $T_\eta$ is a contraction: for any $Z_1, Z_2 \in \mathcal{Z}$, using the nonexpansiveness of Euclidean projections, strong-monotonicity and the Lipschitz property of $F$,

$$ \begin{aligned} \|T_\eta(Z_1) - T_\eta(Z_2)\|^2 &\le \|Z_1 - \eta F(Z_1) - (Z_2 - \eta F(Z_2))\|^2 \\ &= \|Z_1 - Z_2\|^2 - 2\eta \langle Z_1 - Z_2, F(Z_1) - F(Z_2) \rangle + \eta^2 \|F(Z_1) - F(Z_2)\|^2 \le (1 - \mu/L) \|Z_1 - Z_2\|^2. \end{aligned} $$

This implies that the projected forward algorithm with step-size $\eta = 1/L$ is a contractive algorithm. In particular, this includes gradient descent for optimization and gradient-descent-ascent for saddle-point optimization.
## 2.3 Tightness of analysis.

We also show that our analysis of tracking in Theorem 2.1 is tight, in the sense that one cannot significantly improve the tracking bound uniformly for all contractive algorithms. The counter-example is built on a family of quadratic problems and gradient descent. See Appendix B.2 for a proof.

**Theorem 2.2.** For any $C \in (0, 1)$, there exists a sequence of VI problems and a C-contractive algorithm over this sequence such that

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \geq \frac{1}{(1-C)^2} \sum_{t=2}^{T} \|Z_t^* - Z_{t-1}^*\|^2. $$

## 2.4 Tracking guarantees and comparison to existing results.

If a time-varying VI problem is $\alpha$-tame (Definition 1.1), then Theorem 2.1 implies that for a C-contractive algorithm $A$ the tracking error is less than

$$ \tau_T(A) \leq \frac{c}{(1-C)^2} T^{\alpha} + \frac{\|Z_1 - Z_1^*\|^2}{1-C}, $$

where $c > 0$ is the constant in the $\alpha$-tame property.

To our knowledge, our precise setting (exact feedback for time-varying VIs) has not been studied previously at this level of generality. Previous work on time-varying problems has focused either on optimization or games, sometimes in settings with stochastic feedback, in which the learner only gets to observe a noisy value of the operator instead of the true operator. Despite the differences in assumptions, let us provide more detailed points of comparison to the related literature, starting with exact feedback.

Our bound recovers the rates of Mokhtari et al. [2016] when specialized to strongly-convex time-varying optimization. These rates are optimal, in the sense that one cannot improve the dependence on $\alpha$ uniformly over all sequences of strongly-convex functions.

For $\alpha$-tame sequences of games bi-linear zero-sum games, Theorem 8 in Zhang et al. [2022] implies a tracking bound of order $T^{(3+\alpha)/4}$. Our results do not apply in this setting in general, since there are no known contractive algorithms for that class of problems. However, as discussed in Example 2.4, the extra-gradient method is contractive when initialized in a small enough neighborhood of the equilibrium of the game. Hence our results suggest a possible improvement on their bounds for some benign sequences of games, in which the equilibria stay within the basin of attraction of previous equilibria. We leave this investigation for future work. Note also that the algorithm of Zhang et al. [2022] is built on a meta-learning framework, and maintains simultaneously $\log T$ algorithms over $T$ rounds, which may come with an unwieldy memory and computation cost.

Anagnostides et al. [2024] consider time-varying games and minimize the cumulative squared differences between the players' choices and their best responses. Under the additional assumption that each time-varying game is played for multiple rounds, they obtain instance-dependent bounds with respect to the first-order variations of the game plus the variation of the $\epsilon$-Nash equilibria (cf. Theorem 3.5). However, the assumptions restricting the time-variance of the game in this work are not comparable to our setting in Theorem 1.

Finally, tracking guarantees are also studied in the context of time-varying functions or games with stochastic feedback. Note that this setting is fundamentally more challenging and the algorithms and analysis techniques applied in this work are tailored to these specific challenges. In the case of time-varying strongly monotone games, sublinear tracking bounds of order $T^{(2+\alpha)/3}$ are given in Theorem 2 and Corollary 3 in Duvocelle et al. [2023]. Another example for optimization with stochastic feedback is Besbes et al. [2015]: given a tame time-varying strongly convex optimization problem, their Theorem 4 bounds the dynamic regret $\sum_{t=1}^{T} f_t(x_t) - f_t(x_t^*)$ by a term of order $T^{(1+\alpha)/2}$. This directly implies sublinear tracking guarantees. Again, we leave the study of stochastic feedback for future work. The following table compares the most relevant existing results.
**Example 2.2.** Suppose $F: \mathbb{R}^d \to \mathbb{R}^d$ is a maximally $\mu$-strongly monotone operator. Let $R$ denote the resolvent operator, that is, $R = (\mathrm{Id} + F)^{-1}$. (Note that $R$ is related to the backward iteration, and that operator $R$ coincides with the proximal mapping when $F$ is a gradient operator.) Then a point $Z^*$ is a solution to the VIP($F, \mathbb{R}^d$) if and only if $F(Z^*) = 0$, or equivalently if and only if $R(Z^*) = Z^*$. By Theorem 2.1. (xi) in Bauschke et al. [2011] the resolvent operator is $1/(1+\mu)$-contractive This implies in particular that for any $Z \in \mathbb{R}^d$,

$$ \|R(Z) - Z^*\| = \|R(Z) - R(Z^*)\| \le (1 + \mu)^{-1} \|Z - Z^*\|. $$

**Example 2.3.** Consider the following zero-sum game:

$$ \min_{y \in \mathbb{R}} \max_{x \in \mathbb{R}} \left\{ x^2 + 3 \sin^2(x) + a \sin^2(x) \sin^2(y) - y^2 - 3 \sin^2(y) \right\}, \quad (2) $$

where we assume $a \in [0, 1]$. Note that the game in (2) is not a strongly monotone game and it does not satisfy the diagonal concavity assumption of Rosen [1965] (see Proposition B.2 in Appendix B.1 for details). However, note that for any $y \in \mathbb{R}$, the function $\nu(\cdot, y) = x^2 + 3 \sin^2(x) + a \sin^2(x) \sin^2(y) - y^2 - 3 \sin^2(y)$ satisfies the so called *restricted secant inequality* (RSI) (see Proposition B.1, Appendix B.1 for details and definition). We show in Proposition B.3 in Appendix B.1, that contraction is satisfied for gradient descent ascent. A sequence $(a_t)$, with $a_t \in [0, 1]$ defines a time-varying instance of this zero-sum game. Note that this example is inspired and heavily builds on well-known results in optimization by Karimi et al. [2016].

**Example 2.4.** We conclude these examples with a more informal one, regarding the extra-gradient method on bi-affine saddle-point problems, which include in particular equilibria of matrix zero-sum games. Suppose that $\mathcal{Z}$ is the product of two polytopes $\mathcal{P}_1, \mathcal{P}_2$ and that $F$ is a map of the form $F(x, y) = c^\top x + d^\top y - x^\top A y$; that is, $F$ is bi-affine. Korpelevich [1976] introduced the extra-gradient method for saddle-point optimization, which consists in playing

$$ Z_{t+1} = \mathrm{Proj}_{\mathcal{Z}} \{ Z_t - \eta F(\mathrm{Proj}_{\mathcal{Z}} \{ Z_t - \eta F(Z_t) \}) \}. $$

She shows (in Theorem 3 in the reference) that if $F$ has a unique solution $Z_F^*$, then the extra-gradient method is eventually a contraction. Precisely, there exists a (relative) neighborhood of the solution, denoted by $\mathcal{U}(Z_F^*)$, such that extra-gradient becomes a contraction as soon as it reaches $\mathcal{U}(Z_F^*)$. Moreover, suppose in addition the matrix $A$ is a square matrix and the equilibrium lies in the relative interior of both $\mathcal{P}_1$ and $\mathcal{P}_2$, i.e., in $\mathrm{relint}(\mathcal{P}_1) \times \mathrm{relint}(\mathcal{P}_2)$. In that case, the contraction constant depends only on upper bounds of the euclidean operator norms of $A$, and of $A^{-1}$. To fully characterize the contraction property and complete this example, one would need to explicitly describe $\mathcal{U}(Z_F^*)$. Then, given a subset $\mathcal{C} \subset \mathcal{Z}$, the extra-gradient method started at $Z_0 \in \mathcal{C}$ is a contractive algorithm on the set of VI problems VIP($F, \mathcal{Z}$) such that $\mathcal{U}(Z_F^*) \subset \mathcal{C}$. Giving exact descriptions of $\mathcal{U}(Z_F^*)$ is technically challenging and goes beyond the scope of this paper. Recent works (e.g., Wei et al. [2021]) provide finite-time linear convergence bounds with explicit constants, but use arguments that do not explicitly rely upon the contraction property.

## 2.3 Tightness of analysis.

We also show that our analysis of tracking in Theorem 2.1 is tight, in the sense that one cannot significantly improve the tracking bound uniformly for all contractive algorithms. The counter-example is built on a family of quadratic problems and gradient descent. See Appendix B.2 for a proof.

**Theorem 2.2.** For any constant $C \in (0, 1)$, there exists a sequence of VI problems and a $C$-contractive algorithm over this sequence such that

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \ge \frac{1}{(1-C)^2} \sum_{t=2}^{T} \|Z_t^* - Z_{t-1}^*\|^2. $$

## 2.4 Consequences: tracking guarantees.

If a time-varying VI problem is $\alpha$-tame (Definition 1.1), then Theorem 2.1 implies that for a $C$-contractive algorithm $A$ the error is less than

$$ \tau_T(A) \le \frac{c}{(1-C)^2} T^\alpha + \frac{\|Z_1 - Z_1^*\|^2}{1-C}, $$

where $c > 0$ is the constant in the $\alpha$-tame property. In the special case of a time-varying optimization problem, this tracking guarantee together with Section 2.2, Example 2.1 reduce to the results in Mokhtari et al. [2016]. In the context of time-varying zero-sum matrix games, similar tracking guarantees are obtained in Theorem 8 in Zhang et al.
<table><thead><tr><th></th><th>Rate</th><th>Setting</th></tr></thead><tbody><tr><td>Duvocelle et al. [2023]</td><td>T<sup>2+α/3</sup></td><td>Stochastic feedback<br>Strongly monotone time-varying games.</td></tr><tr><td>Besbes et al. [2015]</td><td>T<sup>1+α/2</sup></td><td>Stochastic feedback<br>Strongly convex time-varying functions.</td></tr><tr><td>Zhang et al. [2022]</td><td>T<sup>3+α/4</sup></td><td>Deterministic feedback<br>Time-varying bi-linear games (not strongly monotone).</td></tr><tr><td>Mokhtari et al. [2016]</td><td>T<sup>α</sup></td><td>Deterministic feedback<br>Time-varying strongly convex functions.</td></tr><tr><td>Theorem 2.1</td><td>T<sup>α</sup></td><td>Deterministic feedback<br>Time-varying VIs - any contractive algorithm.</td></tr></tbody></table>

### 3 Periodic time-varying variational inequalities and tracking.

All tracking guarantees of this section are under the assumption that the time-varying variational inequalities have a periodic solution (see Definition 1.2), and that the operators are all $\mu$-strongly monotone. Additional assumptions on $\mathcal{Z}$ or $\mathcal{F}$ are specified at the beginning of each subsection.

**Vacuousness of tracking guarantees via contraction.** We start by showing that the results from Section 2, although applicable, are not sufficient to guarantee tracking for periodic time-varying VI problem. Indeed, consider a $k$-periodic time-varying VI problem. Applying Theorem 2.1, the bound we obtain on $\tau_T(\mathcal{A})$ for $C$-contractive algorithms is

$$ \frac{1}{(1-C)^2} \sum_{t=1}^{T-1} \|Z_t^* - Z_{t+1}^*\|^2 + \frac{\|Z_1 - Z_1^*\|^2}{1-C} \ge \frac{1}{(1-C)^2} \left[ \frac{T}{k} \right]^k \sum_{i=1}^k \|Z_i^* - Z_{i+1}^*\|^2 + \frac{\|Z_1 - Z_1^*\|^2}{1-C} \quad (2) $$

This bound is linear in $T$ as soon $Z_i^* \ne Z_{i+1}^*$ for some $i \in [k]$. However, since the periodic problem is inherently easier than a non-periodic time-varying VI problem, we would expect better tracking guarantees.

Our goal in this section is to improve on the vacuous bound in (2) for periodic time-varying variational inequalities. We propose new algorithms, based on the aggregation of forward-backward methods, that leverage the periodicity of the data. The results of this section are divided into two parts:

1. In Section 3.1, we provide a worst-case regret bound for bounded domains $\mathcal{Z}$, which hold for any sequence of strongly monotone operators $F_t$, against any periodic sequence of comparators (cf. Theorem 3.1). We apply this result in Corollary 3.2 to derive tracking bounds when the comparators are specified to be the solutions of an exact periodic sequence of strongly monotone VI problem.

2. Complementing these results, we show in Section 3.2 that for VIP((F_t), Z) with periodic solutions constant tracking bounds are achievable under the additional assumption that the $F_t$'s are Lipschitz continuous (cf. Theorem 3.2). These results hold even for unbounded domains $\mathcal{Z}$.

#### 3.1 Aggregation of online gradient descent.

We now present our first set of results for tracking solutions of periodic VI problems. At the same time, we introduce the algorithmic framework that will also be used in Section 3.2. In this section, we assume that the domain $\mathcal{Z}$ is bounded, with diameter $D = \max_{Z_1, Z_2 \in \mathcal{Z}} \|Z_1 - Z_2\|$, and that the sequence of operators is bounded by a known constant $G \ge \max_{Z \in \mathcal{Z}} \|F_t(Z)\|$.
[2022]. For α-tame games, that result directly yields a tracking bound of order $T^{(3+\alpha)/4}$. Note that for matrix games, we obtain tracking guarantees due to Example 2.2 Part 4. However, this requires either a fully mixed equilibrium or a warm-up phase to reach the relative neighborhood of the equilibrium. Neither of these assumption is required in Zhang et al. [2022]. However, the algorithm Zhang et al. [2022] apply requires a meta-learner and is therefore computationally more expensive.

Tracking guarantees are also studied in the context of time-varying functions or games with stochastic feedback. Note that due to the stochasticity, this setting is fundamentally more challenging and the algorithms and analysis techniques applied in this work are tailored to these specific challenges. In the case of time-varying strongly monotone games, sublinear tracking bounds of order $T^{2+\alpha/3}$ are obtained in Theorem 2 and Corollary 3 in Duvocelle et al. [2023] allowing for stochastic feedback. Another example of results that allow for stochastic feedback is the work by Besbes et al. [2015]. Given a tame time-varying strongly convex optimization problem, Theorem 4 provides a bound with respect to the dynamic regret $\sum_{t=1}^T f_t(x_t) - f_t(x_t^*)$ of order $T^{(1+\alpha)/2}$. If the time-varying functions are L-smooth, this directly implies sublinear tracking guarantees, that is, sublinear guarantees on $\sum_{t=1}^T \|x_t - x_t^*\|$. It might be an interesting future research direction to investigate if these techniques and our approach could be combined to obtain better results for the more challenging stochastic setting.

# 3 Periodic time-varying variational inequalities and tracking.

The tracking guarantees of this section are under the assumption that the time-varying variational inequalities have a periodic solution (see Definition 1.2), and that the operators are all $\mu$-strongly monotone.

**Vacuousness of tracking guarantees via contraction.** We start by showing that the results from Section 2, although applicable, are not sufficient to guarantee tracking for periodic time-varying VI problem. Indeed, consider a $k$ periodic time-varying VI problem. Applying Theorem 2.1, the bound we obtain on $\tau_T(A)$ for $C$-contractive algorithms is

$$ \frac{1}{(1-C)^2} \sum_{t=1}^{T-1} \|Z_t^* - Z_{t+1}^*\|^2 + \frac{\|Z_1 - Z_1^*\|^2}{1-C} \ge \frac{1}{(1-C)^2} \left[ \frac{T}{k} \right]^k \sum_{i=1}^k \|Z_i^* - Z_{i+1}^*\|^2 + \frac{\|Z_1 - Z_1^*\|^2}{1-C} \quad (3) $$

This bound is linear in $T$ as soon $Z_i^* \ne Z_{i+1}^*$ for some $i \in [k]$. However, since the periodic problem is inherently easier than a non-periodic time-varying VI problem, we would expect better tracking guarantees.

Our goal in this section is to improve on the vacuous bound in (3) for periodic time-varying variational inequalities. We propose new algorithms, based on the aggregation of forward methods, that leverage the periodicity of the data. The results of this section are divided into two parts:

1. In Section 3.1, we provide a worst-case regret bound for bounded domains, which hold for any sequence of strongly monotone operators $F_t$, against any periodic sequence of comparators (cf. Theorem 3.1). We apply this result in Corollary 3.2 to derive tracking bounds when the comparators are specified to be the solutions of an exact periodic sequence of strongly monotone VI problem.

2. Complementing these results, we show in Section 3.2 that for VIP((F_t), Z) with periodic solutions constant tracking bounds are achievable under the additional assumption that the $F_t$'s are Lipschitz continuous (cf. Theorem 3.2). These results hold even for unbounded domains Z.

## 3.1 Aggregation of online gradient descent.

We now present our first set of results for tracking solutions of periodic VI problems. At the same time, we introduce the algorithmic framework that will also be used later in Section 3.2. In this section, we assume that the domain Z is bounded, with diameter $D = \max_{Z_1, Z_2 \in Z} \|Z_1 - Z_2\|$, and that the sequence of operators is bounded by a known constant $G \ge \max_{Z \in Z} \|F_t(Z)\|$.

**On (dynamic) regret bounds for variational inequalities.** The guarantees we derive in this section are upper bounds on quantities of the form

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2, \quad (4) $$

where $\mu \ge 0$ is a positive number, typically the strong-monotonicity constant of the operators. In analogy with the online optimization terminology, we call those guarantees *dynamic regret bounds* for the sequence of operators $(F_t)$,
**On (dynamic) regret bounds for variational inequalities.** The guarantees we derive in this section are upper bounds on quantities of the form

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2, \quad (3) $$

where $\mu \ge 0$ is a positive number, typically the strong-monotonicity constant of the operators. In analogy with the online optimization terminology, we call those guarantees *dynamic regret bounds* for the sequence of operators ($F_t$) against the sequence of comparators ($C_t$). Note that, unless specified otherwise, we make no assumptions on the way the operators $F_t$ are generated: they could be chosen by an adversary. The periodicity of the data is captured by the assumption that the comparators are periodic.

As explained in the proof of Corollary 3.2, when the comparator $C_t$ at time $t$ is specified to be the $k$-periodic solutions to the operators ($F_t$), and when the $F_t$ are $\mu$-strongly monotone, then the dynamic regret is an upper bound on the tracking error.

**Remark 3.1.** Upper bounding the dynamic regret (3) can be useful beyond the application we give. For example, in the two fundamental examples of convex optimization and convex-concave saddle-point optimization, under an additional strong monotonicity assumption, the dynamic regret (3) is an upper bound on the natural measure of performance, as explained in Examples 3.1 and 3.2.

**Remark 3.2.** Reducing the computation of solutions of VI problems to the optimization of scalar objectives is standard in the VI literature. In fact, the summand in (3) is a special case of a saddle-function associated to the operator $F_t$, as discussed extensively in Larsson and Patriksson [1994]; precisely, the saddle-function $\mathcal{L}: (z, c) \mapsto \langle F(z), z - c \rangle - (\mu/2)\|z - c\|^2$. This saddle-function naturally defines a 'gap function', $\sup_{c \in \mathcal{X}} \mathcal{L}(z, c)$, see also Facchinei and Kanzow [2007] (Definition 10.2.2), which can serve as an optimization proxy for the VI.

**Example 3.1 (Optimization).** If $f$ is a $\mu$-strongly-convex function, with gradient $F = \nabla f$, then by strong convexity, for any $z, c \in \mathcal{Z}$,

$$ f(z) - f(c) \le \langle F(z), z - c \rangle - \frac{\mu}{2} \|z - c\|^2. $$

**Example 3.2 (Saddle-point Optimization).** If $u: (p_1, p_2) \mapsto u(p_1, p_2)$ is a $\mu$-strongly-convex-strongly-concave function, with pseudo-gradient $F = (\nabla_1 u, -\nabla_2 u)$, then for any $z = (z_1, z_2)$ and $c = (c_1, c_2) \in \mathcal{Z}$, we have:

$$ \begin{aligned} u(z_1, c_2) - u(c_1, z_2) &= u(z_1, c_2) - u(z_1, z_2) + u(z_1, z_2) - u(c_1, z_2) \\ &\le -\langle \nabla_2 u(z_1, z_2), z_2 - c_2 \rangle - \frac{\mu}{2} \|z_2 - c_2\|^2 + \langle \nabla_1 u(z_1, z_2), z_1 - c_1 \rangle - \frac{\mu}{2} \|z_1 - c_1\|^2 \\ &= \langle F(z), z - c \rangle - \frac{\mu}{2} \|z - c\|^2. \end{aligned} $$

In these two central examples, regret guarantees imply upper bounds on the cumulative gaps against the sequence of comparators ($C_t$)$_{t \in \mathbb{N}}$. Note also that this notion has previously appeared in the online learning literature, e.g. in Zhang et al. [2022].

**General structure of aggregation and base algorithm.** We use a *meta-algorithm template*, as stated in Algorithm 2. The basic idea is to maintain several copies of a *base algorithm* with different parameters. The outputs of the base algorithms are combined using an aggregation algorithm, guaranteeing performance comparable to that of the best base algorithm. We refer to this aggregation algorithm as the *meta algorithm*. In our setting, the base algorithms are initialized with different period lengths $i \in [K]$, i.e., $A_{cFB}^{(1)}, \dots, A_{cFB}^{(K)}$, and the meta algorithm guarantees a performance comparable to an algorithm initialized with the correct cycle length $k$, $A_{cFB}^{(k)}$. This idea is commonly used in online learning for adapting to unknown parameters in step-size tuning, e.g., Van Erven et al. [2021]. The unknown parameter in our case is the period $k$ of the sequence.

**Base algorithm.** The base algorithm we consider is an instance of gradient descent tailored to periodic problems, instead of stationary problems. Given a period $i$, the algorithm maintains $i$ independent iterates and updates them cyclically. We denote by $\tilde{F}_t$ the operator given as input to the algorithm: this is to emphasize the fact that we will sometimes give a surrogate operator to the algorithm instead of the true operator $F_t$ (see (4)); we also instantiate it with the true operator in Section 3.2.

We first recall the analysis for a fixed comparator. This is a direct application of the online gradient descent bound on the sequence of linear losses $z \mapsto \langle \tilde{F}_t(Z_t), z \rangle$, used in the analysis of online gradient descent for strongly-convex losses; see Theorem 2.1 in Bartlett et al. [2007].
against the sequence of comparators ($C_t$). Note that, unless specified otherwise, we make no assumptions on the way the operators $F_t$ are generated: they could be chosen by an adversary. The periodicity of the data is captured by the assumption that the comparators are periodic. Such regret bounds are of interest for multiple reasons.

First, as explained in the proof of Corollary 3.2, when the comparator $C_t$ at time $t$ is specified to be the $k$-periodic solutions to the operators ($F_t$), and when the $F_t$ are $\mu$-strongly monotone, the dynamic regret is an upper bound on the tracking error.

Moreover, variational inequalities often come equipped with a natural measure of the merit of proposed solutions. Bach and Levy [2019] identify settings in which this measure has nice properties and call these settings variational inequalities with a gap function; the gap function ($x, y) \mapsto \Delta_F(x, y)$ quantifies the quality of $x$ compared to $y$ as a solution of the VIP($F, \mathcal{Z}$). We recall the definition for completeness and provide examples showing that (4) is often a natural upper bound on the gap function. Although we technically never use the concept of gap functions in our results, we add this small detour to highlight the relation to existing literature.

**Definition 3.1 (Bach and Levy [2019]).** Given a variational inequality problem ($F, \mathcal{Z}$), we say the function $\Delta_F: \mathcal{Z} \times \mathcal{Z}$ is a gap function for $F$ if it is convex in its first argument, if

$$ \forall Z_1, Z_2 \in \mathcal{Z}: \quad \Delta_F(Z_1, Z_2) \le \langle F(Z_1), Z_1 - Z_2 \rangle, $$

and if $Z^*$ is a solution of $F$ if and only if $\Delta_F(Z^*, Z) \le 0$ for all $Z \in \mathcal{Z}$.

Two fundamental examples are convex optimization and convex-concave saddle-point optimization. In both cases, under an additional strong monotonicity assumption, the regret (4) for the VI is an upper bound on the gap function between $z_t$ and $c_t$. Let us expand on these examples.

**Example 3.1 (Optimization).** If $f$ is a $\mu$-strongly-convex function, with gradient $F = \nabla f$, then the quantity $\Delta_F(x, y) = f(x) - f(y)$ is a gap function for $F$. Moreover, by strong convexity, for any $z, c \in \mathcal{Z}$, we have

$$ \Delta_F(z, c) = f(z) - f(c) \le \langle F(z), z - c \rangle - \frac{\mu}{2} \|z - c\|^2. $$

**Example 3.2 (Saddle-point Optimization).** If $u: (p_1, p_2) \mapsto u(p_1, p_2)$ is a $\mu$-strongly-convex-strongly-concave function, with pseudo-gradient $F = (\nabla_1 u, -\nabla_2 u)$, then the quantity $\Delta_F((p_1, p_2), (q_1, q_2)) = u(p_1, q_2) - u(q_1, p_2)$ is a gap function for $F$, and for any $z = (z_1, z_2)$ and $c = (c_1, c_2) \in \mathcal{Z}$, we have:

$$ \begin{aligned} \Delta_F(z, c) &= u(z_1, c_2) - u(c_1, z_2) = u(z_1, c_2) - u(z_1, z_2) + u(z_1, z_2) + u(c_1, z_2) \\ &\le -\langle \nabla_2 u(z_1, z_2), z_2 - c_2 \rangle - \frac{\mu}{2} \|z_2 - c_2\|^2 + \langle \nabla_1 u(z_1, z_2), z_1 - c_1 \rangle - \frac{\mu}{2} \|z_1 - c_1\|^2 \\ &= \langle F(z), z - c \rangle - \frac{\mu}{2} \|z - c\|^2. \end{aligned} $$

In these two central examples, regret guarantees imply upper bounds on the cumulative gaps against the sequence of comparators ($c_t$)$_{t \in \mathbb{N}}$. Note also that this notion has previously appeared in the online learning literature, e.g. in Zhang et al. [2022].

**General structure of aggregation and base algorithm.** We use a meta-algorithm template, as stated in Algorithm 2. The basic idea is to maintain several copies of a base algorithm with different parameters. The outputs of the base algorithms are combined using an aggregation algorithm, guaranteeing performance comparable to that of the best base algorithm. This idea is commonly used in online learning for adapting to unknown parameters in step-size tuning, e.g., Van Erven et al. [2021]. The unknown parameter in our case is the period $k$ of the sequence.

**Base algorithms.** The base algorithm we consider is an instance of gradient descent tailored to periodic problems, instead of stationary problems. Given a period $i$, the algorithm maintains $i$ independent iterates and updates them cyclically. We denote by $\tilde{F}_t$ the operator given as input to the algorithm: this is to emphasize the fact that we will sometimes give a surrogate operator to the algorithm instead of the true operator $F_t$ (see (5)); we also instantiate it with the true operator in Section 3.2.

We first recall the analysis for a fixed comparator. This is a direct application of the online gradient descent bound on the sequence of linear losses $z \mapsto \langle \tilde{F}_t(Z_t), z \rangle$, used in the analysis of online gradient descent for strongly-convex losses; see Theorem 2.1 in Bartlett et al. [2007].

**Proposition 3.1.** Given $\mu > 0$, for any fixed $C \in \mathcal{Z}$, if $\eta_t = 1/(\mu t)$, then the Online Forward Method (Algorithm 1 tuned with period $i = 1$), ensures that for all $t \in [T]$,

$$ \sum_{t=1}^{T} \langle \tilde{F}_t(Z_t), Z_t - C \rangle - \frac{\mu}{2} \|Z_t - C\|^2 \le \frac{(\max_{t \in [T]} \|\tilde{F}_t(Z_t)\|)^2}{2\mu} (\log(T) + 1). $$
**Algorithm 1:** Algorithm $A_{cFB}^{(i)}$: Cyclic Forward-Backward Method

**Input:** Assumed period $i$, step-size schedule ($\eta_s$), flag if used as a base algorithm.

**Initialization:** Initial iterates ($Z_{1,1}, Z_{2,1}, \dots, Z_{i,1}$) = ($Z_1, Z_1, \dots, Z_1$) with $Z_1 \in \mathbb{Z}$.

**for** $t = 1, \dots, T$ **do**

Pick current index $n = (t \mod i) + 1$, and corresponding time $s = \lceil t/i \rceil$;

Play $Z_t = Z_{n,s}$, receive operator $\tilde{F}_t$;

if $A_{cFB}^{(i)}$ functions as a base algorithm then

Forward $Z_t$ to meta-algorithm;

Update $Z_{n,s+1} = \text{Proj}_{\mathbb{Z}}(Z_{n,s} - \eta_s \tilde{F}_t(Z_{n,s}))$;

**Output:** sequence ($Z_t$);

**Proposition 3.1.** Given $\mu > 0$, for any fixed $C \in \mathbb{Z}$, if $\eta_t = 1/(\mu t)$, then the Online Forward-Backward Method (Algorithm ?? tuned with period $i = 1$), ensures that for all $t \in [T]$,

$$ \sum_{t=1}^{T} \langle \tilde{F}_t(Z_t), Z_t - C \rangle - \frac{\mu}{2} \|Z_t - C\|^2 \le \frac{(\max_{t \in [T]} \|\tilde{F}_t(Z_t)\|)^2}{2\mu} (\log(T) + 1). $$

A direct consequence of this bound is that if we tune the learner with a period $k$, then for any $k$-periodic sequence of comparators ($C_t$)$_{t \in \mathbb{N}}$ with $C_t = C_{t+k}$ for all $t \in \mathbb{N}$, the regret of cyclic gradient descent can be upper bounded by separating the update times of each individually maintained iterate.

**Corollary 3.1.** Let ($F_t$)$_{t \in \mathbb{N}}$ be a sequence of G-bounded operators. For any $k$-periodic sequence of comparators ($C_t$)$_{t \in \mathbb{N}}$ with $C_t = C_{t+k}$ for all $t \in \mathbb{N}$, Algorithm ?? tuned with the period $k$ and step-size schedule $\eta_s = 1/(\mu s)$ played on $\tilde{F}_t = F_t$ enjoys the bound

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 \le \frac{kG^2}{2\mu} (\log(T/k) + 1). $$

We note that this result holds in particular if the sequence of comparators is the (k-periodic) sequence of solutions, that is $C_t = Z_t^*$. Since ($F_t$)$_{t \in \mathbb{N}}$ are by assumption $\mu$-strongly monotone, this implies that the tracking error is upper bounded by the regret. Thus, $\tau_T(A_{cFB}^{(k)})$ is in the order of $(G^2/\mu^2)k \log T$. In practice, we expect the period of the data to be unknown. Therefore, we run several instances of the base algorithm. Each instance is initialized with a different period length, and the iterates of the base algorithms are aggregated by a meta-algorithm.

**Aggregation.** We maintain $K \ge k$ base algorithms and aggregate them using an expert algorithm. We wish to retain a regret bound for the combined aggregated iterates that is comparable to that of the best algorithm in hindsight. The standard Hedge bound is of order $\sqrt{T \ln K}$, which would give an overhead cost that dominates the bound of the best-performing algorithm. In order to retain the fast logarithmic rate, we propose a carefully crafted aggregation scheme working for VI problems, based on a combination of surrogate operators with exponentially weighted aggregation (a.k.a. Vovk's Aggregating Algorithm in this context, see Vovk [1995]). The $K$ base algorithms are updated with the affine surrogate operators

$$ \tilde{F}_t : z \mapsto F_t(Z_t) + \mu(z - Z_t). \quad (4) $$

The term 'surrogate operators' is inspired by the online learning terminology. From each base algorithm $A_{cFB}^{(i)}$, the meta-algorithm receives an iterate $Z_t^{(i)}$. Then, the meta-algorithm combines the base predictions $Z_t^{(i)}$ by maintaining a probability distribution over the $K$ base algorithms. This probability distribution is given by the aggregation algorithm. We define for all $i \in [K]$,

$$ p_{t,i} \propto \exp\left(-\lambda \sum_{s=1}^{t-1} \ell_{s,i}\right) \quad \text{with} \quad \ell_{s,i} = \langle F_s(Z_s), Z_s^{(i)} \rangle + (\mu/2) \|Z_s^{(i)} - Z_s\|^2. \quad (5) $$

The learning rate is set to $\lambda = (4\mu(D + G/\mu)^2)^{-1}$. This learning rate depends on the diameter of the set $\mathbb{Z}$, i.e., $D = \max_{Z,Z' \in \mathbb{Z}} \|Z - Z'\|$, and $G \ge \max_{t \in T} \|F_t(Z_t)\|$ and the strong-monotonicity parameter $\mu$. We then select the final action $Z_t = p_{t,1}Z_t^{(1)} + \dots + p_{t,K}Z_t^{(K)}$.
**Algorithm 1:** Algorithm A<sub>cFW</sub><sup>(i)</sup>: Cyclic Forward Method

**Input:** Assumed period *i*, step-size schedule ($\eta_s$).

**Initialization:** Initial iterates ($Z_{1,1}$, $Z_{2,1}$, ..., $Z_{i,1}$) = ($Z_1$, $Z_1$, ..., $Z_1$) with $Z_1 \in \mathbb{Z}$.

**for** *t* = 1, ..., *T* **do**

Pick current index *n* = (*t* mod *i*) + 1, and corresponding time *s* = ⌊*t*/*i*⌋;

Play $Z_t = Z_{n,s}$, receive operator $\tilde{F}_t$;

Update $Z_{n,s+1} = \text{Proj}_{\mathbb{Z}}(Z_{n,s} - \eta_s \tilde{F}_t(Z_{n,s}))$;

A direct consequence of this bound is that if we tune the learner with a period *k*, then for any *k*-periodic sequence of comparators ($\tilde{C}_t$)$_{t \in [k]}$, the regret of cyclic gradient descent can be upper bounded by separating the update times of each individually maintained iterate.

**Corollary 3.1.** For any *k*-periodic sequence of comparators ($C_t$), Algorithm 1 tuned with the period *k* and step-size schedule $\eta_s = 1/(\mu s)$ played on $\tilde{F}_t = F_t$ enjoys the bound

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 \le \frac{kG^2}{2\mu} (\log(T/k) + 1). $$

We note that this result holds in particular if the sequence of comparators is the (k-periodic) sequence of solutions, that is $C_t = Z_t^*$. Since ($F_t$)$_t$ are by assumption $\mu$-strongly monotone, this implies that the tracking error is upper bounded by the regret. Thus, $\tau_T(A_{cFW}^{(k)})$ is in the order of $k \log T$, up to a multiplicative factor of $G^2/\mu^2$. In practice, we expect the period of the data to be unknown. Therefore, we run several instances of the base algorithm. Each instance is initialized with a different period length and the iterates of the base algorithms are aggregated by a meta-algorithm.

**Aggregation.** We maintain $K \ge k$ base algorithms and aggregate them using an expert algorithm. We wish to retain a regret bound for the combined aggregated iterates that is comparable to that of the best algorithm in hindsight. The standard Hedge bound is of order $\sqrt{T \ln K}$, which would give an overhead cost that dominates the bound of the best-performing algorithm. In order to retain the fast logarithmic rate, we propose a carefully crafted aggregation scheme working for VI problems, based on a combination of surrogate operators with exponentially weighted aggregation (a.k.a. Vovk's Aggregating Algorithm in this context, see Vovk [1995]). The K base algorithms are updated with the affine surrogate operators

$$ \tilde{F}_t : z \mapsto F_t(Z_t) + \mu(z - Z_t). \qquad (5) $$

From each base algorithm $A_{cFW}^{(i)}$, the meta-algorithm receives an iterate $Z_t^{(i)}$. Then, the meta-algorithm combines the base predictions $Z_t^{(i)}$ by maintaining a probability distribution over the K base algorithms. This probability distribution is given by the aggregation algorithm. We define for all $i \in [K]$,

$$ p_{t,i} \propto \exp\left(-\lambda \sum_{s=1}^{t-1} \ell_{s,i}\right) \quad \text{with} \quad \ell_{s,i} = \langle F_s(Z_s), Z_s^{(i)} \rangle + (\mu/2) \|Z_s^{(i)} - Z_s\|^2. \qquad (6) $$

The learning rate is set to $\lambda = (4\mu(D + G/\mu)^2)^{-1}$. This learning rate depends on the diameter of the set $\mathbb{Z}$, i.e., $D = \max_{Z,Z' \in \mathbb{Z}} \|Z - Z'\|$, and $G \ge \max_{t \in T} \|F_t(Z_t)\|$ and the strong-monotonicity parameter $\mu$. We then select the final action $Z_t = p_{t,1}Z_t^{(1)} + \cdots + p_{t,K}Z_t^{(K)}$.

**Algorithm 2:** Meta-Algorithm Template

**Data:** Maximum cycle period $K \in \mathbb{N}$, aggregator over K base iterates, K base algorithms.

**for** *t* = 1, ..., *T* **do**

Receive base plays ($Z_t^{(1)}$, ..., $Z_t^{(K)}$) from base algorithms;

Receive distribution $p_t$ over $[K]$ base algorithms from aggregator algorithm;

Play $Z_t = p_{t,1}Z_t^{(1)} + \cdots + p_{t,K}Z_t^{(K)}$ and receive $F_t$;

Update aggregator algorithm with losses $\ell_{i,t}$;

Update base algorithms with operator $\tilde{F}_t$;
**Algorithm 2:** Meta-Algorithm Template

**Data:** Maximum cycle period $K \in \mathbb{N}$, aggregator over $K$ base iterates, $K$ base algorithms.
**for** $t = 1, \dots, T$ **do**

Receive base plays $(Z_t^{(1)}, \dots, Z_t^{(K)})$ from base algorithms (Algorithm ??);

Receive distribution $p_t$ over $[K]$ base algorithms from aggregator algorithm;

Play $Z_t = p_{t,1}Z_t^{(1)} + \dots + p_{t,K}Z_t^{(K)}$ and receive $F_t$;

Update aggregator algorithm with losses $\ell_{i,t}$;

Update base algorithms (Algorithm ??) with operator $\tilde{F}_t$;

**Theorem 3.1.** Let $(F_t)_{t \in \mathbb{N}}$ be a sequence of $G$-bounded operators on the bounded domain $\mathcal{Z}$ of diameter $D$. Consider Algorithm 2, initialized with maximum period $K$, using exponential weights (with losses specified in (5), and learning rate $\lambda = (4\mu(D+G/\mu)^2)^{-1}$) and cyclic forward-backward method (Algorithm ??) with step-size schedule $\eta_s = 1/(\mu s)$ as base algorithms. For any $k \le K$, against any $k$-periodic sequence of comparators $(C_t)_{t \in \mathbb{N}}$ with $C_t = C_{t+k}$ for all $t \in \mathbb{N}$

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 \le \frac{(G + \mu D)^2}{2\mu} \left( k \log(T/k) + k + 8 \log K \right). $$

Interestingly, since the mappings $z \mapsto \langle F_t(z), z \rangle$ are not convex, our proof requires careful manipulations of the surrogates. A positive computational side-effect of using these surrogates is that the full procedure only requires a single evaluation of the operator $F_t$. For a complete proof, see Appendix C.1.

**Corollary 3.2 (Tracking periodic variational inequalities).** Let $(F_t)_{t \in \mathbb{N}}$ be a sequence of $\mu$-strongly monotone and $G$-bounded operators on the bounded domain $\mathcal{Z}$ of diameter $D$. If the solutions of $(F_t)$ are $k$-periodic, then Algorithm 2 tuned as in the statement of Theorem 3.1 guarantees

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le \frac{(G + \mu D)^2}{\mu^2} \left( k \log(T/k) + k + 8 \log K \right). $$

*Proof.* Proof. Note that $\langle F(Z_t^*), Z_t - Z_t^* \rangle \ge 0$ by definition. With $\mu$-strong monotonicity, we have

$$ \mu \|Z_t - Z_t^*\|^2 \le \langle F_t(Z_t) - F(Z_t^*), Z_t - Z_t^* \rangle \le \langle F_t(Z_t), Z_t - Z_t^* \rangle, $$

which we rearrange to make the instantaneous regret appear

$$ \frac{\mu}{2} \|Z_t - Z_t^*\|^2 \le \langle F_t(Z_t), Z_t - Z_t^* \rangle - \frac{\mu}{2} \|Z_t - Z_t^*\|^2. $$

Therefore, after summing over $t \in [T]$, Theorem 3.1 provides the claimed result.

## 3.2 Constant tracking error under Lipschitzness.

The previous results are limited to the constrained setting, in which the domain $\mathcal{Z}$ is bounded. They also require the knowledge of an upper bound $G$ on the norm of the point values of the operators. We now consider the unconstrained setting, $\mathcal{Z} = \mathbb{R}^d$, and add a Lipschitz assumption on the operators $F_t$. Note that we still assume that the operators are all $\mu$-strongly monotone. Furthermore, in contrast to Section 3.1, we allow for $K$ evaluations of the function $F_t$. Our algorithm is another instance of the meta-algorithm framework introduced in Section 3.1. In this Section, we give an algorithm that obtains a constant tracking bound that is independent of $T$.

**Base algorithms.** To take advantage of the Lipschitz assumption, we use variants of the cyclic forward algorithm with a constant learning rate of $\eta = 1/L$ and feedback $\tilde{F}_t = F_t$. Note that using this feedback now requires $K$ evaluations of $F_t$. However, the constant learning rate together with full multi-point feedback gives exponentially fast convergence to the solution.

**Lemma 3.1.** Let VIP((F_t), Z) be a time-varying VI problem with $k$-periodic solutions. Assume the functions $F_t$ are $\mu$-strongly monotone and L-Lipschitz. Let $(Z_t^{(k)})$ denote the sequence of iterates of the forward algorithm $A_{c\text{FB}}^{(k)}$ (Algorithm ??), tuned with $\eta_t = 1/L$ and correct period $k$. Then

$$ \|Z_t^{(k)} - Z_t^*\|^2 \le \left(1 - \frac{\mu}{L}\right)^{[t/k]} \|Z_1 - Z_t^*\|^2. $$
**Theorem 3.1.** Consider Algorithm 2, tuned with maximum period K, using exponential weights (with losses specified in (6), and learning rate $\lambda = (4\mu(D + G/\mu)^2)^{-1}$) and cyclic forward method (Algorithm 1) with step-size schedule $\eta_s = 1/(\mu s)$ as base algorithms. For any $k \le K$, against any k-periodic sequence of comparators $(C_t)$,

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 \le \frac{(G + \mu D)^2}{2\mu} \left( k \log(T/k) + k + 8 \log K \right). $$

Interestingly, since the mappings $z \mapsto \langle F_t(z), z \rangle$ are not convex, our proof requires careful manipulations of the surrogates. A positive computational side-effect of using these surrogates is that the full procedure only requires a single evaluation of the operator $F_t$. For a complete proof, see Appendix C.1.

**Corollary 3.2 (Tracking periodic variational inequalities).** If the sequence of $\mu$-strongly monotone operators $(F_t)_{t \in \mathbb{N}}$ has a k-periodic solution for some $k \le K$, then Algorithm 2 tuned as in the statement of Theorem 3.1 guarantees

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le \frac{(G + \mu D)^2}{\mu^2} \left( k \log(T/k) + k + 8 \log K \right). $$

*Proof.* Note that $\langle F(Z_t^*), Z_t - Z_t^* \rangle \ge 0$ by definition. With $\mu$-strongly monotonicity, we have

$$ \mu \|Z_t - Z_t^*\|^2 \le \langle F_t(Z_t) - F(Z_t^*), Z_t - Z_t^* \rangle \le \langle F_t(Z_t), Z_t - Z_t^* \rangle, $$

which we rearrange to make the instantaneous regret appear

$$ \frac{\mu}{2} \|Z_t - Z_t^*\|^2 \le \langle F_t(Z_t), Z_t - Z_t^* \rangle - \frac{\mu}{2} \|Z_t - Z_t^*\|^2. $$

Therefore, after summing over $t \in [T]$, Theorem 3.1 provides the claimed result.

<figure><img src="image_3.png" /></figure>

## 3.2 Constant tracking error under Lipschitzness.

The previous results are limited to the constrained setting, in which the domain $\mathcal{Z}$ is bounded. They also require the knowledge of an upper bound $G$ on the norm of the point values of the operators. We now consider the unconstrained setting, $\mathcal{Z} = \mathbb{R}^d$, and add a Lipschitz assumption on the operators $F_t$. Note that we still assume that the operators are all $\mu$-strongly monotone. Furthermore, in contrast to Section 3.1, we allow for $K$ evaluations of the function $F_t$. Our algorithm is another instance of the meta-algorithm framework introduced in Section 3.1. In this Section, we give an algorithm that obtains a constant tracking bound that is independent of $T$.

**Base algorithms.** To take advantage of the Lipschitz assumption, we use variants of the cyclic forward algorithm with a constant learning rate of $\eta = 1/L$ and feedback $\tilde{F}_t = F_t$. Note that using this feedback now requires $K$ evaluations of $F_t$. However, the constant learning rate together with full multi-point feedback gives exponentially fast convergence to the solution.

**Lemma 3.1.** Let VIP((F_t), Z) be a time-varying VI problem with k-periodic solutions. Assume the functions F_t are $\mu$-strongly monotone and L-Lipschitz. Let $(Z_t^{(k)})$ denote the sequence of iterates of the forward algorithm $A_{\text{CFW}}^{(k)}$ (Algorithm 1), tuned with $\eta_t = 1/L$ and correct period k. Then

$$ \|Z_t^{(k)} - Z_t^*\|^2 \le \left(1 - \frac{\mu}{L}\right)^{\lfloor t/k \rfloor} \|Z_1 - Z_t^*\|^2. $$

The result is a small modification of Example 2.1 in Section 2.2. For details, see Appendix C.2.1. From this, we obtain that the base algorithm tuned with the correct period $k$ converges exponentially fast to the solutions $Z_1^*, \dots, Z_k^*$.

**Aggregation with adaptive step-sizes.** The distribution over the base algorithms is an instance of exponential weights, with the crucial difference with respect to our approach of Section 3.1 that we use an adaptive learning rate $\lambda_t$. More precisely,

$$ p_{t,i} \propto \exp\left(-\lambda_t \sum_{s=1}^{t-1} \ell_{s,i}\right) \quad \text{with} \quad \ell_{s,i} = \langle F_s(Z_s), Z_s^{(i)} \rangle + (\mu/2) \|Z_s^{(i)} - Z_s\|^2. \quad (7) $$
The result is a small modification of Example 2.1 in Section 2.2. From this, we obtain that the base algorithm tuned with the correct period k converges exponentially fast to the solutions $Z_1^*, \dots, Z_k^*$.

**Aggregation with adaptive step-sizes.** The distribution over the base algorithms is an instance of exponential weights, with the crucial difference with respect to our approach of Section 3.1 that we use an adaptive learning rate $\lambda_t$. More precisely, let $p_1$ be the uniform distribution over $[K]$. Define the loss vector and the loss of the average play at time $t$,

$$ \ell_{t,i} = \langle F_t(Z_t), Z_t^{(i)} \rangle + (\mu/2) \|Z_t^{(i)} - Z_t\|^2 \quad \text{and} \quad \bar{\ell}_t = \langle F_t(Z_t), Z_t \rangle. \qquad (6) $$

For the first rounds, if $\bar{\ell}_t \le \min_i \ell_{t,i}$, define $p_{t+1}$ to be the uniform distribution over the minimal components of $\ell_1 + \dots + \ell_t$. Let $T_0$ be the first time at which $m_t < \bar{\ell}_t$. For $t \ge T_0 + 1$, denoting by $(x)_+ = \max(x, 0)$ the positive part of $x$, set $p_t$ to be the distribution over $[K]$ such that

$$ p_{t,i} \propto \exp\left(-\lambda_t \sum_{s=1}^{t-1} \ell_{s,i}\right) \quad \text{where} \quad \lambda_t = \frac{\log K}{\sum_{s=T_0}^{t-1} (\bar{\ell}_s - m_s)_+} \quad \text{and} \quad m_s = -\frac{1}{\lambda_s} \sum_{i=1}^{K} p_{t,i} e^{-\lambda_s \ell_{s,i}}. \qquad (7) $$

Note that the denominator in the definition of the learning rate is positive, by definition of $T_0$. The initialization steps, when $t \le T_0$ can be informally described as setting $\lambda_t = \infty$ when that denominator is 0. The learning rate is tuned in the style of AdaHedge (De Rooij et al. [2014]), to decrease proportionally to the inverse of the cumulative gap between the loss of the average of plays $\bar{\ell}_s$, and the mix loss $m_s$.

In the unconstrained setting, the final convergence guarantees unavoidably depend on the initial point chosen by the algorithm. Given a time-varying problem VIP((F_t), Z) with periodic solutions, let $\mathcal{C}$ denote the set of all the solutions of the operators (F_t) (which has cardinality at most k). For simplicity, we assume that all base algorithms are initialized with the same $Z_1 \in \mathcal{Z}$, and define the initial distance

$$ D_0 = \max_{Z^* \in \mathcal{C}} \|Z_1 - Z^*\|. $$

The next theorem provides constant tracking guarantees for our algorithm in the unconstrained setting, depending on the initial distance and condition number $\kappa = L/\mu$.

**Theorem 3.2.** Consider a time-varying unconstrained problem VIP((F_t), $\mathbb{R}^d$) with k-periodic solutions. Assume all operators (F_t) are all $\mu$-strongly monotone and L-Lipschitz, and let $\kappa = L/\mu$ denote the condition number of the problem. Algorithm 2, tuned with $K \ge k$ and losses and learning rate specified in (6) and (7), has tracking guarantees:

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le 4D_0^2(2+\kappa)(2(2\kappa^2+1)(2+\kappa)\log K + (2\kappa+1)\kappa k + 1). $$

The proof is detailed in Appendix C.2. The crux of the proof, in Lemmas C.1 and C.2, consists in showing that thanks to the adaptive learning rate, the aggregation method ensures similar regret bounds as in the proof of Theorem 3.1, but replacing $D$ and $G$ by their empirical counterparts $\hat{D}_T$ and $\hat{G}_T$, defined as follows:

$$ \hat{D}_T = \max_{t \in [T], i,j \in [K]} \left\| Z_t^{(i)} - Z_t^{(j)} \right\| \quad \text{and} \quad \hat{G}_T = \max_{t \in [T], i \in [K]} \left\| F_t(Z_t^{(i)}) \right\|. $$

we show in Lemma C.3 in Appendix C.2 that these quantities are controlled by $D_0$ for our algorithms and assumptions. The rest of the proof is similar to that of Theorem 3.1.

# 4 Application Examples

In this section, we provide examples of time-varying VI problems to illustrate the practical relevance of our results. As mentioned in the introduction, it is well known that VIs are a generalization of optimization and saddle-point problems. Examples 4.1 and 4.2 in Section 4.1 are illustrations of practical problems in these settings. See also Facchinei and Pang [2003, Chapter I] for many more examples, which often have natural time-varying extensions. In Section 4.2, we give an example illustrating that VIs can also be applied to statistical parameter estimation. In its full generality, this problem is not equivalent to an optimization problem, thus showcasing the wider scope of VIs.
We tune the learning rate in the style of AdaHedge (De Rooij et al. [2014]), to decrease proportionally to the inverse of the cumulative gap between the loss of the average of play $\bar{l}_s$, and the mix loss $m_s$, given below. Denoting by $(x)_+ = \max(x, 0)$ the positive part of a real number $x$, define

$$ \lambda_t = \frac{\log K}{\sum_{s=1}^{t-1} (\bar{l}_s - m_s)_+}, \quad \text{where} \quad \bar{l}_s = \langle F_s(Z_s), Z_s \rangle \quad \text{and} \quad m_s = -\frac{1}{\lambda_s} \sum_{i=1}^{K} p_{t,i} e^{-\lambda_s l_{s,i}}. \qquad (8) $$

When $\sum_{s=1}^{t-1} (\bar{l}_s - m_s)_+ = 0$, we use the natural conventions that $\lambda_t = \infty$, that $p_t$ is the uniform distribution over the minimal components of the vector $\ell_1 + \cdots + \ell_{t-1}$, and that $m_t = \min_{i \in [K]} \ell_{t,i}$.

In the unconstrained setting, the final convergence guarantees unavoidably depend on the initial point chosen by the algorithm. Given a time-varying problem VIP((F_t), Z) with periodic solutions, let $\mathcal{C}$ denote the set of all the solutions of the operators (F_t) (which has cardinality at most k), and define the initial distance. For simplicity, we assume that all base algorithms are initialized with the same $Z_1 \in Z$.

$$ D_0 = \max_{Z^* \in \mathcal{C}} \|Z_1 - Z^*\|. $$

The next theorem provides constant tracking guarantees for our algorithm in the unconstrained setting, depending on the initial distance and condition number $\kappa = L/\mu$.

**Theorem 3.2.** Consider a time-varying unconstrained problem VIP((F_t), R^d) with k-periodic solutions. Assume all operators (F_t) are all $\mu$-strongly monotone and L-Lipschitz, and let $\kappa = L/\mu$ denote the condition number of the problem. Algorithm 2, tuned with $K \ge k$ and losses and learning rate specified in (7) and (8), has tracking guarantees:

$$ \sum_{t=1}^{T} \|Z_t - Z_t^*\|^2 \le 4D_0^2(2 + \kappa) (2(2\kappa^2 + 1)(2 + \kappa) \log K + (2\kappa + 1)\kappa k + 1). $$

The proof is detailed in Appendix C.2. The crux of the proof, in Lemmas C.1 and C.2, consists in showing that thanks to the adaptive learning rate, the aggregation method ensures similar regret bounds as in the proof of Theorem 3.1, but replacing D and G by empirical counterparts $\hat{D}_T$ and $\hat{G}_T$, defined as follows:

$$ \hat{D}_T = \max_{t \in [T], i, j \in [K]} \|Z_t^{(i)} - Z_t^{(j)}\| \quad \text{and} \quad \hat{G}_T = \max_{t \in [T], i, j \in [K]} \|F_t(Z_t^{(i)})\|. $$

we show in Lemma C.3 in Appendix C.2 that these quantities are controlled by $D_0$ for our algorithms and assumptions. The rest of the proof is similar to that of Theorem 3.1.

**Limitations and future work.** The algorithm of Section 3.1 is computationally more efficient than the algorithm in Section 3.2. However, in contrast to the latter, it requires a bounded domain and provides (only) logarithmic guarantees. It would be an interesting question for future work to combine the merits of both approaches. Another limitation becomes clear when considering the special case of games: In this case, the aggregation schemes discussed in this section may be unnatural, as they require coordination between the players. Decoupling the aggregation in general cases would be an interesting direction for further investigation.

We also remark that both results in Section 3.1 and 3.2 require deterministic feedback. An interesting direction is to consider more robust models, in which we would observe the inexact values of the operators. Examples of such models could be a stochastic feedback model or time-varying VI problems that admit periodic behavior only in the limit. These feedback models would capture application examples, such as fluctuating markets, in a more realistic setting. Inexact feedback models are inherently more challenging, the convergence guarantees in this setting are weaker, and the aggregation step would be more difficult to implement.

# 4 Main negative result: chaos.

In this section, we provide a detailed case study of the convergence properties of fixed step-size gradient descent on some periodic sequences problems. We focus on optimization instead of the general VI framework since we aim to provide simple examples of striking phenomena. For the same reason, we focus on periodic functions. Recall that this implies that the functions have a periodic solution (cf. Definition 1.2). Accordingly, we always go for the smallest dimension in which these observations occur.

Besides the surprising appearance of chaotic behavior which we found interesting in its own, we draw two practical conclusions from the results in this section:
## 4.1 Examples: Time-varying convex optimization and saddle-point problems.

For clarification and completeness, we recall the relation between VIs and convex optimization and saddle-point problems. Let $\mathcal{Z}$ be a closed convex set and let $f: \mathcal{Z} \to \mathbb{R}$ be a differentiable and convex function. In the following examples, we will either assume that $\mathcal{Z}$ is bounded or that $f$ is strongly convex. Under either of these assumptions, if we set $F = \nabla f$ then $u^*$ is a solution of VIP($F, \mathcal{Z}$) if and only if

$$u^* \in \operatorname*{argmin}_{x \in \mathcal{Z}} f(x). \qquad (\text{OPT})$$

If the function $f$ is $\mu$-strongly convex and $L$-smooth, then the corresponding VI problem is $\mu$-strongly monotone and $L$-Lipschitz. To illustrate the relation to game theory, consider an $n$-player normal-form game with convex loss functions: Let $\mathcal{Z}^{(i)}$ denote the (convex) action set of player $i \in [n]$, and define $\mathcal{Z} = \prod_{i \in [n]} \mathcal{Z}^{(i)}$, and $\mathcal{Z}^{(-i)} = \prod_{i \in [n] \setminus \{i\}} \mathcal{Z}^{(i)}$. Finally, denote by $\nu^{(i)}: \mathcal{Z}^{(i)} \times \mathcal{Z}^{(-i)} \to \mathbb{R}$ the $i$-th player's loss function and assume that for all $x^{(-i)} \in \mathcal{Z}^{(-i)}$, the loss $\nu^{(i)}(\cdot, x^{(-i)})$ is convex and differentiable. A point $\bar{x} = [\bar{x}^{(1)}, \dots, \bar{x}^{(n)}] \in \mathcal{Z}$ is a Nash equilibrium if

$$\forall i \in [n], \forall x^{(i)} \in \mathcal{Z}^{(i)}: \quad \nu^{(i)}(x^{(i)}, \bar{x}^{(-i)}) \geq \nu^{(i)}(\bar{x}^{(i)}, \bar{x}^{(-i)}). \qquad (\text{NE})$$

Consider the operator $F$ mapping $(x^{(1)}, \dots, x^{(n)})$ to $[\nabla_{x^{(1)}} \nu^{(1)}(x^{(1)}, x^{(-1)}), \dots, \nabla_{x^{(n)}} \nu^{(n)}(x^{(n)}, x^{(-n)})]$, also called the pseudo-gradient of the game $(\nu^{(i)})$. The set of Nash equilibria is exactly the set of solutions of VIP($F, \mathcal{Z}$). See, e.g., Rockafellar and Wets [1998, Chapter 12], or Kinderlehrer and Stampacchia [2000, Chapter 7], for more details and examples.

**Example 4.1 (Repeated kelly auctions on seasonal markets).** We adapt an example from Duvocelle et al. [2023]. Consider a repeated Kelly auction (Kelly et al. [1998]) where $n$ bidders participate in a repeated auction on a limited good. Let $R_t \ge 0$ denote the minimal price, i.e., the entry barrier set at time $t$. For all $i \in [n]$, given a maximum budget $b_t^{(i)}$, and their value for the good, $v_t^{(i)} \ge 0$, player $i$ bids $x_t^{(i)} \in [0, b_t^{(i)}]$. Each player gets a fraction of the good relative to their bid. Their utility function is

$$\nu_t^{(i)}(x_t^{(i)}, x_t^{(-i)}) = v_t^{(i)} \frac{x_t^{(i)}}{R_t + \sum_{j \in [n]} x_t^{(j)}} - x_t^{(i)}.$$

For many goods, the demand or the supply depends on exogenous factors that vary with time. For example, the availability of agricultural products depends on the season, and the demand for luxury goods may depend on fashion trends. These conditions naturally affect the minimal price $R_t$, and the individual values $v_t^{(i)}$, making the game time-varying with seasonal shifts. For further examples and context, also see Section 2 in Rahme et al. [2021].

To apply the results from Section 2, assume that the demand adheres to a slowly shifting fashion trend, thus satisfying tameness (cf. Definition 1.1). Furthermore, note that the pseudo-gradient is Lipschitz and monotone; thus, adding a quadratic regularizer to the utility guarantees strong monotonicity. Hence, combining Example 2.1 and Theorem 2.1 provides tracking guarantees for this example. Further, consider a scenario where the availability of a good follows an idealized seasonal trend, hence satisfying periodicity. Again, consider the regularized utilities. Given that we have an upper estimation on the periodicity, Theorem 3.2 provides a constant tracking bound. We note that the assumption of periodicity is rather unrealistic for many real-world applications. Hence, this example also illustrates the practical relevance of the future research directions in Section 3.2.

**Example 4.2 (Linear regression on streaming data).** Consider a linear regression problem with a continuous stream of i.i.d. data as introduced in Foster [1991]. That is, let $\{(a_1, b_1), \dots, (a_{n_t}, b_{n_t})\} \subset \mathbb{R}^d \times \mathbb{R}$ denote the data at time $t$, and $A_t := [a_1, \dots, a_{n_t}]^\top \in \mathbb{R}^{n_t \times d}$ and $b_t := [b_1, \dots, b_{n_t}]^\top \in \mathbb{R}^{n_t}$ the corresponding data matrix and vector, respectively. Suppose the data points $(a_i, b_i)$ are coming from a stream of measurements of the same experiment. In an experimental setting, it can be beneficial to have a real-time regression on the data stream (see e.g. the software tool [Mat, 2020]). Hence, define $f_t: \mathbb{R}^d \to \mathbb{R}$ with $f_t(x) = \|A_t x - b_t\|^2 + \lambda \|x\|^2$. Under the assumption that the experimental data follows a fixed sub-gaussian distribution, the solution to the empirical regression problem approaches the solution with respect to the true distribution with high probability, see Shalev-Shwartz and Ben-David [2014, Section 13]. This defines an optimization problem with a slowly shifting solution and thus a *tame time-varying* optimization problem.

We note that the objective function of this tame time-varying optimization problem is $\lambda$-strongly convex; hence, combining Example 2.2 and Theorem 2.1 provides tracking guarantees for this example. Note that Example 2.2 does not require $\mathcal{Z}$ to be bounded. Hence, it fits this example where $x \in \mathbb{R}^d$ is from an unbounded set.
* The tuning of the parameters for the forward algorithm (i.e., gradient descent) is crucial: while the forward algorithm for fixed strongly monotone and Lipschitz continuous VI can guarantee exponential convergence with a fixed step size, the combination of ignoring the periodicity and step size tuning may lead to undesirable behavior in the periodic case.

* In finite-sum optimization problems, when using gradient descent and going through multiple passes on the functions, common wisdom regarding the tuning of the learning rate does not apply. Indeed, in some cases, even if the iterates are in an unstable regime, increasing the learning rate may boost performance (see the example in the introduction and Observation 1.)

## 4.1 Notation

Let $x_1 \in \mathbb{Z}$ be an arbitrary point and consider the sequence of iterates generated by gradient descent with step-size $\eta$ on a $k$-periodic sequence of $\mu$-strongly convex and $L$-smooth problems $(f_t)_{t \in \mathbb{N}}$, and define $F_t = \nabla f_t$. Define the Gradient Descent iteration operator

$$ \Phi_{\eta,t} : x \mapsto \mathrm{Proj}_{\mathbb{Z}}(x - \eta F_t(x)) \quad \text{and} \quad \bar{\Phi}_{\eta} = \Phi_{\eta,1} \circ \cdots \circ \Phi_{\eta,k}. \qquad (9) $$

With this notation, the sequence $(x_t)_{t \in \mathbb{N}}$ satisfies the equations $x_{t+1} = \Phi_{\eta,t}(x_t)$ for all $t \in \mathbb{N}$, and $x_{(n+1)k} = \bar{\Phi}_{\eta}(x_{nk})$ for all $k \in \mathbb{N}$.

Since all $\Phi_{\eta,t}$ are continuous (by continuity of Euclidean projections on convex sets, together with the $L$-Lipschitzness of $F_t$), and since $\mathbb{Z}$ is closed, the sequence $(x_t)_{t \in \mathbb{N}}$ is bounded if and only if one of its subsequences $(x_{i+kn})_n$ is bounded. Therefore, the boundedness of the iterates is solely dictated by $\bar{\Phi}_{\eta}$, and the starting $F_1$ does not matter for this discussion. We already provided an example of quadratic time-varying problems with a constant minimizer—arguably the simplest setting. We observed that the value of the step-size $\eta$ with the best speed of convergence may be larger than values that lead to divergence. In the following section, we make a series of striking observations, notably:

* The sequence of GD iterates can converge to arbitrary points, go through cycles of arbitrary length, or even be chaotic, see Section 4.3.

* Period-halving bifurcations occur, interleaved with divergent phases, see Figure 1

* For specific choices of the step size, we observe star-shaped limit sets, see Section 4.4.

## 4.2 Examples of exotic behavior with large step-sizes.

We describe simple instances of time-varying problems which exhibit exotic behavior of gradient descent with large step sizes. We show that even in simple examples the trajectory of Gradient Descent —and its dependence on the learning rate— can be surprisingly complex. In particular, we show that chaos can emerge in one-dimensional periodically time-varying optimization problems with a fixed minimizer. Moreover, increasing the learning rate from a chaotic trajectory may lead to divergence, to cyclic behavior, or to converging trajectories.

The examples we consider are built using the following instance of a time-varying optimization problem. Given a periodic sequence of symmetric positive definite matrices $A_t$, let $f_t$ be functions from $\mathbb{R}^d$ to $\mathbb{R}$ defined by

$$ f_t(x) = \log(1 + e^{x^\top A_t x/2}), \quad \text{so that} \quad F_t(x) := f'_t(x) = \frac{e^{x^\top A_t x/2}}{1 + e^{x^\top A_t x/2}} A_t x. \qquad (10) $$

The functions $f_t$ are close to quadratic functions (and $F_T$ are close to linear), but the non-linearity in the gradients is enough to produce chaotic behavior. Note also that the minimizer is always 0. Furthermore, note that $(F_t)$ is a sequence of $\mu$-strongly monotone and $L$-Lipschitz functions, where the parameters $\mu$ and $L$ depend on $(A_t)$.

**Definitions: stability and chaos.** Our description of the trajectory of gradient descent uses notions from real discrete dynamical systems theory which we now recall. We refer the reader to the two monographs Block and Coppel [2006], Ruette [2017] from which we extracted these definitions and theorems. Note that we do not state the results in full generality, for the sake of brevity.

Following the notations introduced above, we denote a trajectory of gradient descent with step-size $\eta$ by $(x_t)$ (the step-size will be clear from context). That is, $x_{t+1} = x_t - \eta F_t(x_t)$, and $x_{k(n+1)} = \bar{\Phi}_{\eta}(x_{kn})$ and the points at time steps that are multiples of the time-period $k$ are orbits of the stationary dynamical system $(\bar{\Phi}_{\eta}, \mathbb{R})$. (Note that $\bar{\Phi}_{\eta}$ is always continuous.)

A point $x$ is said to be periodic under the map $\Phi$ if its (forward) orbit $(\Phi^n(x))_n$ is periodic; the period of $x$ is the smallest period of the sequence. A periodic orbit $P = \{z_1, \dots, z_p\}$ of $\Phi$ is asymptotically stable if there exists an
## 4.2 Example: Generalized Linear Models

To illustrate the generality of VIs and their time-varying variant, we conclude this section with an example from statistical parameter estimation.

**Example 4.3 (Generalized linear models).** We present an example from Juditsky and Nemirovsky [2019]. Consider a generalized linear model corrupted by Gaussian noise; i.e., we obtain observations $(a_t, b_t) \in \mathbb{R}^d \times \mathbb{R}$, $t = 1, \dots, T$, that satisfy $b_t = \phi(\langle Z^*, a_t \rangle) + \xi_t$, where $\phi: \mathbb{R} \mapsto \mathbb{R}$ is a continuous and increasing link function, $\xi_t \sim \mathcal{N}(0, \sigma^2)$, and $Z^* \in \mathbb{R}^d$ is an unknown parameter vector. Juditsky and Nemirovsky [2019] observe that this parameter estimation problem can be reformulated as a monotone VI problem VIP($F, \mathbb{R}^d$) with operator

$$F(Z) = \frac{1}{T} \sum_{t=1}^{T} a_t [\phi(\langle Z, a_t \rangle) - b_t]. \qquad (8)$$

Moreover, under minimal distributional assumptions, it can be proven that the population-type operator (i.e., the expected value of $F$ under the target distribution) becomes strongly monotone. For simplicity of presentation, we focus here on the univariate case where the link function $\phi$ maps from $\mathbb{R} \to \mathbb{R}$. In this specific setting, the operator $F$ in (8) is always the gradient of a convex potential function, making the VI equivalent to a convex optimization problem. However, Juditsky and Nemirovsky [2019] consider a broader setting with multivariate link functions $(\phi: \mathbb{R}^m \to \mathbb{R}^m)$, where the operator's Jacobian is not necessarily symmetric. Consequently, the general problem is not equivalent to an optimization problem and serves as a genuine example of the broader VI framework. The above VI can also be considered in a streaming setting, such as proposed in Example 4.2. If needed, regularization can be introduced by adding $\lambda Z$ to the operator in equation (8).

Note that similarly to Example 4.2, Theorem 2.1 combined with Example 2.2 imply tracking guarantees.

**Remark 4.1.** The above VI reformulation is one common option for solving GLM. Conversely, we may also consider the maximum likelihood estimator of the unknown parameter $Z \in \mathbb{R}^d$ via solving the (nonlinear) least squares problem

$$\min_{Z \in \mathbb{R}^d} \left\{ \frac{1}{T} \sum_{t=1}^{T} [b_t - \phi(\langle Z, a_t \rangle)]^2 \right\}.$$

This optimization problem is typically nonconvex (with the exception of the case where the link function is affine). Thus, computing a solution might be computationally expensive.

# 5 Main negative result: chaos.

In this section, we provide a detailed case study of the convergence properties of fixed step-size gradient descent on some periodic time-varying problems. We focus on optimization instead of the general VI framework since we aim to provide simple examples of striking phenomena. For the same reason, we focus on periodic functions. Recall that this implies that the functions have a periodic solution (cf. Definition 1.2). Accordingly, we always aim for the smallest dimension in which these observations occur.

**Motivation and practical relevance** Besides the surprising appearance of chaotic behavior, which we found interesting in its own, we draw several practical conclusions from the results in this section. First, proper parameter tuning is crucial: while fixed step sizes ensure exponential convergence in strongly monotone and Lipschitz continuous variational inequalities, ignoring periodicity in step size tuning can lead to instability. Second, in many applications, optimization methods are used with step sizes exceeding the bounds for which theoretical convergence results exist. Empirical studies (see Cohen et al. [2021]) show that larger step sizes can improve machine learning model performance in both convergence and generalization. Our study adds a surprising piece to this puzzle: in periodic settings, increasing the step size can turn a divergent system into a convergent one. Given the popularity of single-shuffle stochastic gradient descent [Mishchenko et al., 2020], which falls into the periodic setting, this finding might have practical implications. However, we emphasize that the emergence of chaos, as discussed here, is a property of the autonomous (time-invariant) dynamical system induced by the periodic composition of operators. In general non-autonomous (arbitrary time-varying) or stochastic settings (using SGD), Li-Yorke chaos is not meaningful. Moreover, in practice, adaptive step-size methods (such as AdaGrad, AcceleGrad, UnixGrad, UnderGrad, Levy et al. [2018], Kavis et al. [2019], Antonakopoulos et al. [2022] or Armijo-type line-search Armijo [1966]) are often employed and can mitigate instability or chaotic behavior in deterministic settings. Our results thus primarily highlight a theoretical phenomenon in deterministic, periodic settings, and serve as a caution for using fixed step-sizes in such cases. Investigating whether similar phenomena occur more generally in stochastic systems, or with adaptive step sizes is left for future work.
open set $U$ such that for any $x \in U$, the sequence of repeated iterations $(\Phi^{(pn)}(x))_n$ converges to a point in $P$. Then any trajectory initialized in that open set will eventually get close to the orbit $P$.

The results we are most interested in apply to cases in which the function $\Phi$ maps a compact interval $I \subset \mathbb{R}$ into itself. Then, we shall say that the pair $(\Phi, I)$ forms an interval map. In that case, the orbits from a starting point inside $I$ stay bounded in $I$.

Given an interval map $\Phi$ on $I$, a set $S \subset I$ is said to be *scrambled* if for any $x, y \in S$, we have $\limsup_{n \to \infty} |\Phi^n(x) - \Phi^n(y)| > 0$ and $\liminf_{n \to \infty} |\Phi^n(x) - \Phi^n(y)| = 0$. The interval map $(\Phi, I)$ is said to be *chaotic* in the sense of Li-Yorke if there exists an uncountable scrambled set.

Li-Yorke chaos is one of the standard formalizations of chaos in real-valued dynamical systems, and has been invoked recently in the game theory and optimization literature to describe the trajectories of standard algorithms in games, e.g. in Bielawski et al. [2021].

## 4.3 Different behaviours of gradient descent.

In our first example, we consider $d = 1$ and a periodic sequence of length 2. More precisely, in (10) let for all $n \in \mathbb{N}$,

$$A_{2n+1} = 0.25 \quad \text{and} \quad A_{2n} = 4.$$

Despite the very basic features of the problem, we show that depending on the value of the learning rate, Gradient Descent can:

- converge to zero; the global minimizer of the sequence,

- oscillate around cycles of points which are not minimizers,

- be chaotic in the sense of Li-Yorke,

- diverge to infinity.

Surprisingly, increasing the learning rate can get the trajectory out of unstable dynamics and recover convergence. This is summarized in the bifurcation diagram of Figure 1, and in Observation 1 below.

<figure><img src="image_4.png" /><figcaption>Figure 1: Bifurcation diagram for the dynamical system generated by gradient descent on the periodic problem (10), with the learning rate as the varying parameter, and with initial point $x = -0.1$.</figcaption></figure>

The diagram represents the accumulation points of the sequence of iterates of gradient descent. Grey areas correspond to values of the learning rate for which the sequence diverged. Details on the simulation that lead to this figure can be found in Appendix D.

While an exhaustive theoretical analysis justifying the bifurcation diagram is out of reach (even in the well-studied case of the logistic map), we discuss some specific values of $\eta$ in the following result. There, we put together a collection of
## 5.1 Notation

Let $x_1 \in \mathbb{Z}$ be an arbitrary point and consider the sequence of iterates generated by gradient descent with step-size $\eta$ on a $k$-periodic sequence of $\mu$-strongly convex and $L$-smooth problems ($f_t)_{t \in \mathbb{N}}$, and define $F_t = \nabla f_t$. Define the Gradient Descent iteration operator

$$ \Phi_{\eta,t} : x \mapsto \mathrm{Proj}_{\mathbb{Z}}(x - \eta F_t(x)) \quad \text{and} \quad \bar{\Phi}_{\eta} = \Phi_{\eta,1} \circ \cdots \circ \Phi_{\eta,k}. \qquad (9) $$

With this notation, the sequence $(x_t)_{t \in \mathbb{N}}$ satisfies the equations $x_{t+1} = \Phi_{\eta,t}(x_t)$ for all $t \in \mathbb{N}$, and $x_{(n+1)k} = \bar{\Phi}_{\eta}(x_{nk})$ for all $k \in \mathbb{N}$.

Since all $\Phi_{\eta,t}$ are continuous (by continuity of Euclidean projections on convex sets, together with the $L$-Lipschitzness of $F_t$), and since $\mathbb{Z}$ is closed, the sequence $(x_t)_{t \in \mathbb{N}}$ is bounded if and only if one of its subsequences $(x_{i+kn})_n$ is bounded. Therefore, the boundedness of the iterates is solely dictated by $\bar{\Phi}_{\eta}$, and the starting $F_1$ does not matter for this discussion. We already provided an example of quadratic time-varying problems with a constant minimizer—arguably the simplest setting. We observed that the value of the step-size $\eta$ with the best speed of convergence may be larger than values that lead to divergence. In the following section, we make a series of striking observations, notably:

* The sequence of GD iterates can converge to arbitrary points, go through cycles of arbitrary length, or even be chaotic, see Section 5.3.

* Period-halving bifurcations occur, interleaved with divergent phases, see Figure 1

* For specific choices of the step size, we observe star-shaped limit sets, see Section 5.4.

## 5.2 Examples of chaotic behavior with large step-sizes.

We describe simple instances of time-varying problems which exhibit exotic behavior of gradient descent with large step sizes. We show that even in simple examples the trajectory of Gradient Descent —and its dependence on the learning rate— can be surprisingly complex. In particular, we show that chaos can emerge in one-dimensional periodically time-varying optimization problems with a fixed minimizer. Moreover, increasing the learning rate from a chaotic trajectory may lead to divergence, to cyclic behavior, or to converging trajectories.

The examples we consider are built using the following instance of a time-varying optimization problem. Given a periodic sequence of symmetric positive definite matrices $A_t$, let $f_t$ be functions from $\mathbb{R}^d$ to $\mathbb{R}$ defined by

$$ f_t(x) = \log (1 + e^{x^\top A_t x/2}), \quad \text{so that} \quad F_t(x) := f'_t(x) = \frac{e^{x^\top A_t x/2}}{1 + e^{x^\top A_t x/2}} A_t x. \qquad (10) $$

The functions $f_t$ are close to quadratic functions (and $F_T$ are close to linear), but the non-linearity in the gradients is enough to produce chaotic behavior. Note also that the minimizer is always 0. Furthermore, note that $(F_t)$ is a sequence of $\mu$-strongly monotone and $L$-Lipschitz functions, where the parameters $\mu$ and $L$ depend on $(A_t)$.

**Definitions: stability and chaos.** Our description of the trajectory of gradient descent uses notions from real discrete dynamical systems theory which we now recall. We refer the reader to the two monographs Block and Coppel [2006], Ruette [2017] from which we extracted these definitions and theorems. Note that we do not state the results in full generality, for the sake of brevity.

Following the notations introduced above, we denote a trajectory of gradient descent with step-size $\eta$ by $(x_t)$ (the step-size will be clear from context). That is, $x_{t+1} = x_t - \eta F_t(x_t)$, and $x_{k(n+1)} = \bar{\Phi}_{\eta}(x_{kn})$ and the points at time steps that are multiples of the time-period $k$ are orbits of the stationary dynamical system $(\bar{\Phi}_{\eta}, \mathbb{R})$. (Note that $\bar{\Phi}_{\eta}$ is always continuous.)

A point $x$ is said to be periodic under the map $\Phi$ if its (forward) orbit $(\Phi^n(x))_n$ is periodic; the period of $x$ is the smallest period of the sequence. A periodic orbit $P = \{z_1, \dots, z_p\}$ of $\Phi$ is asymptotically stable if there exists an open set $U$ such that for any $x \in U$, the sequence of repeated iterations $(\Phi^{(pn)}(x))_n$ converges to a point in $P$. Then any trajectory initialized in that open set will eventually get close to the orbit $P$.

The results we are most interested in apply to cases in which the function $\Phi$ maps a compact interval $I \subset \mathbb{R}$ into itself. Then, we shall say that the pair $(\Phi, I)$ forms an interval map. In that case, the orbits from a starting point inside $I$ stay bounded in $I$.

Given an interval map $\Phi$ on $I$, a set $S \subset I$ is said to be *scrambled* if for any $x, y \in S$, we have $\limsup_{n \to \infty} |\Phi^n(x) - \Phi^n(y)| > 0$ and $\liminf_{n \to \infty} |\Phi^n(x) - \Phi^n(y)| = 0$. The interval map $(\Phi, I)$ is said to be *chaotic* in the sense of Li-Yorke if there exists an uncountable scrambled set.
numerical observations and proofs, providing detailed evidence that gradient descent can behave in arbitrarily complex ways on simple time-varying optimization problems.

**Observation 1.** Consider the trajectories of gradient descent with fixed step size applied to the 2-periodic sequence of functions (10). Depending on the value of $\eta$, the trajectories may have any of the following behaviors:

- For $\eta \le 1/2$ or $\eta \in [7.6, 8.4]$, there exists a neighbourhood $\mathcal{U}$ of 0 such that GD converges to 0 for any initial point $x_0 \in \mathcal{U}$; moreover $-0.1 \in \mathcal{U}$.

- For $\eta = 2$ the sequence diverges to $\infty$ for any non-zero initial point.

- For $\eta = 3.9$, there exists an 8-periodic sequence ($c_n$) and an open set $\mathcal{U}$ containing 0.1, such that for any $x_0 \in \mathcal{U}$, we have $(x_n - c_n) \to 0$.

- For $\eta = 6.1$, the system is Li-Yorke chaotic. Precisely, there exists two uncountable sets $S_e$ and $S_o$ such that $S_e$ is scrambled by $\bar{\Phi}_\eta$ and $S_o$ is scrambled by $\Phi_{\eta,2} \circ \Phi_{\eta,1}$.

We substantiate those observation by combining numerical evidence from simulations with classical theorems from discrete dynamical systems theory. In particular, we exhibit a period-3 orbit of $\bar{\Phi}_\eta$, which implies the chaotic behavior by the Period-Three theorem of Li and Yorke [1975], by examining the functions $I - \bar{\Phi}_\eta$ and $I - (\bar{\Phi}_\eta)^3$.

*Proof.* Evidence and proofs substantiating Observation 1.

$\eta \le 1/2$ or $\eta \in [7.6, 8.4]$. The second derivative of the map $f: x \mapsto \log(1 + e^{x^2/2})$ is the function $(1 + e^{-x^2/2} + x^2e^{-x^2/2})/(1 + e^{-x^2/2})^2$, which takes values in $[1/2, 1.31)$. Therefore $f_2(x) = f(2x)$ is 2-strongly-convex and $(4 \times 1.31)$-smooth; similarly, $f_1(x) = f(x/2)$ is $(1/8)$-strongly-convex and $(1.31/4)$-smooth.

We recall that for any $\bar{\mu}$-strongly convex, and $\bar{L}$-smooth function minimized at 0, we have for any $x \in \mathbb{R}^d$, for any $\alpha \ge 0$,

$$ \begin{aligned} \|x - \alpha \nabla f(x)\|^2 &= \|x\|^2 + \alpha^2 \|\nabla f(x)\|^2 - 2\alpha \langle x, \nabla f(x) \rangle \\ &= \|x\|^2 + \alpha^2 \|\nabla f(x) - \nabla f(0)\|^2 - 2\alpha \langle x - 0, \nabla f(x) - \nabla f(0) \rangle \\ &\le \|x\|^2 + \bar{L}^2 \|x\|^2 - 2\bar{\mu} \|x\|^2 = (1 - 2\bar{\mu}\alpha + \alpha^2 \bar{L}^2) \|x\|^2. \end{aligned} $$

In particular for any $n \in \mathbb{N}$,

$$ \|x_{2n+3}\|^2 \le (1 - 2\eta\mu/4 + \eta^2 L/4) \|x_{2n+2}\|^2 \le (1 - 2\eta(4\mu) + \eta^2(4L))(1 - 2\eta\mu/4 + \eta^2 L/4) \|x_{2n+1}\|^2, $$

where $\mu = 1/2$ and $L \le 1.31$. The value of the factor is less than 0.3 for the step-size $\eta = 1/2$. Therefore the sequence converges to 0, regardless of the initial point.

Note moreover that 0 is always a fixed point of $\bar{\Phi}_\eta$, and that $|\bar{\Phi}_\eta)'(0)| < 1$ for any $\eta \in [7.6, 8.4]$. In that case, 0 is an attracting fixed point, proving the claim. Simulations show that $-0.1$ belongs to the basin of attraction of that point.

$\eta = 2$. For $\eta = 2$, we have $|\bar{\Phi}_\eta(x)| > 2|x|$ for any $x \in \mathbb{R} \setminus \{0\}$ (see Figure 2); therefore the trajectory diverges exponentially fast to $\infty$ from any initial point.

$\eta = 3.9$. The map $\psi = (\bar{\Phi}_\eta)^{(4)}$ admits multiple fixed points; see Figure 3(a). For example, Newton's method for computing fixed points (i.e., the sequence defined by $x_{n+1} = x_n - \psi(x_n)/\psi'(x_n)$), initialized at $x = -0.1$, converges numerically to a point $x_\infty \approx -1.35$; its orbit under $\bar{\Phi}_\eta$ is the cycle with approximate values $\{-1.35, 5.92, -1.57, 7.04\}$. The product of $(\bar{\Phi}_\eta)'(u)$ over $u$ in this cycle is approximately $-0.26$, with absolute value less than 1. Therefore, by the next result, this cycle is asymptotically stable for $\bar{\Phi}_\eta$.

**Theorem 4.1** (Prop. 20 Chap. 5 in Block and Coppel [2006]). Let $P = \{x_1, \dots, x_p\}$ be a periodic orbit of a map $\Phi$, and assume that $\Phi$ is differentiable at all $x_i \in P$. If $|\Phi'(x_1) \dots \Phi'(x_p)| < 1$, then $P$ is asymptotically stable.

This shows that there exists an asymptotically stable 4-cycle for $\bar{\Phi}_\eta$, i.e., that the sequence $(x_{2n+1})$ will converge to that cycle when initialized in the appropriate (open) set. The result for the whole sequence follows by noting that $x_{2n+2} = \Phi_{\eta,1}(x_{2n+1})$ for any $n \in \mathbb{N}$, and that $\Phi_{\eta,1}$ is continuous.
Li-Yorke chaos is one of the standard formalizations of chaos in real-valued dynamical systems, and has been invoked recently in the game theory and optimization literature to describe the trajectories of standard algorithms in games, e.g. in Bielawski et al. [2021].

## 5.3 Different behaviours of gradient descent.

In our first example, we consider $d = 1$ and a periodic sequence of length 2. More precisely, in (10) let for all $n \in \mathbb{N}$,

$$A_{2n+1} = 0.25 \quad \text{and} \quad A_{2n} = 4.$$

Despite the very basic features of the problem, we show that depending on the value of the learning rate, Gradient Descent can:

- converge to zero; the global minimizer of the sequence,

- oscillate around cycles of points which are not minimizers,

- be chaotic in the sense of Li-Yorke,

- diverge to infinity.

Surprisingly, increasing the learning rate can get the trajectory out of unstable dynamics and recover convergence. This is summarized in the bifurcation diagram of Figure 1, and in Observation 1 below.

<figure><img src="image_3.png" /><figcaption>Figure 1: Bifurcation diagram for the dynamical system generated by gradient descent on the periodic problem (10), with the learning rate as the varying parameter, and with initial point $x = -0.1$.</figcaption></figure>

The diagram represents the accumulation points of the sequence of iterates of gradient descent. Grey areas correspond to values of the learning rate for which the sequence diverged, axis were broken to hide large areas of divergence. Details on the simulation used to create this figure can be found in Appendix E.

While an exhaustive theoretical analysis justifying the bifurcation diagram is out of reach (even in the well-studied case of the logistic map), we discuss some specific values of $\eta$ in the following result. There, we put together a collection of numerical observations and proofs, providing detailed evidence that gradient descent can behave in arbitrarily complex ways on simple time-varying optimization problems.

**Observation 1.** Consider the trajectories of gradient descent with fixed step size applied to the 2-periodic sequence of functions (10). Depending on the value of $\eta$, the trajectories may have any of the following behaviors:
<figure><img src="image_5.png" /><figcaption>Figure 2: Evidence that $\bar{\Phi}_{\eta}(x) > 2|x|$ for all $x \in \mathbb{R}$, when $\eta = 2$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 3: Evidence of the existence of periodic points for $\bar{\Phi}_{\eta}$. Each intersection of the plotted curves with the x-axis corresponds to a fixed point of the iterated function, and thus of a periodic orbit for $\bar{\Phi}_{\eta}$. A fixed point of $\Phi_{\eta}^{(3)}$ that is not a fixed point of $\Phi_{\eta}$ has period 3; the existence of such points implies chaos.</figcaption></figure>

$\eta = 6.1$. The map $\bar{\Phi}_{\eta}$ is an interval map on $[-2.5, 2.5]$. There exists a 3-periodic point $x \approx 0.20$; its orbit is approximately $\{-0.10, 0.20, 0.04\}$. (See Figure 3(b).) Therefore, by the Period-Three theorem [Li and Yorke, 1975], there exists an uncountable scrambled set $S$ in $[-2.5, 2.5]$ for $\bar{\Phi}_{\eta}$. This corresponds to the odd times of the gradient descent trajectory. The result for even times follows again by considering the fact $x_{2n+2} = \Phi_{\eta,1}(x_{2n+1})$, and that $\Phi_{\eta,1}$ is continuous.

<figure><img src="image_7.png" /></figure>

## 4.4 Iterated function systems and star-shaped attractors.

Another interesting phenomenon we observe for time-varying periodic functions is the convergence to a star-shaped set.

**Example 4.1** (Periodic time-varying problems with star-shaped limit sets). We consider again a 2-periodic sequence as defined in (10), here with dimension $d = 2$ and

$$A_1 = \begin{bmatrix} 3/4 & 0 \\ 0 & 5 \end{bmatrix}, \quad A_2 = \begin{bmatrix} 5 & 1 \\ 1 & 3/4 \end{bmatrix}.$$
- For $\eta \le 1/2$ or $\eta \in [7.6, 8.4]$, there exists a neighbourhood $\mathcal{U}$ of 0 such that GD converges to 0 for any initial point $x_0 \in \mathcal{U}$; moreover $-0.1 \in \mathcal{U}$.

- For $\eta = 2$ the sequence diverges to $\infty$ for any non-zero initial point.

- For $\eta = 3.9$, there exists an 8-periodic sequence $(c_n)$ and an open set $\mathcal{U}$ containing 0.1, such that for any $x_0 \in \mathcal{U}$, we have $(x_n - c_n) \to 0$.

- For $\eta = 6.1$, the system is Li-Yorke chaotic. Precisely, there exists two uncountable sets $S_e$ and $S_o$ such that $S_e$ is scrambled by $\bar{\Phi}_\eta$ and $S_o$ is scrambled by $\Phi_{\eta,2} \circ \Phi_{\eta,1}$.

We substantiate those observations by combining numerical evidence from simulations with classical theorems from discrete dynamical systems theory. In particular, we exhibit a period-3 orbit of $\bar{\Phi}_\eta$, which implies the chaotic behavior by the Period-Three theorem of Li and Yorke [1975], by examining the functions $I - \bar{\Phi}_\eta$ and $I - (\bar{\Phi}_\eta)^3$.

*Proof.* Evidence and proofs substantiating Observation 1.

$\eta \le 1/2$ or $\eta \in [7.6, 8.4]$. The second derivative of the map $f: x \mapsto \log(1 + e^{x^2/2})$ is the function $(1 + e^{-x^2/2} + x^2e^{-x^2/2})/(1 + e^{-x^2/2})^2$, which takes values in $[1/2, 1.31)$. Therefore $f_2(x) = f(2x)$ is 2-strongly-convex and $(4 \times 1.31)$-smooth; similarly, $f_1(x) = f(x/2)$ is $(1/8)$-strongly-convex and $(1.31/4)$-smooth. Therefore, using the same derivation as in Example 2.1, for any $n \in \mathbb{N}$,

$$ \|x_{2n+3}\|^2 \le (1 - 2\eta\mu/4 + \eta^2 L/4)\|x_{2n+2}\|^2 \le (1 - 2\eta(4\mu) + \eta^2(4L))(1 - 2\eta\mu/4 + \eta^2 L/4)\|x_{2n+1}\|^2, $$

where $\mu = 1/2$ and $L \le 1.31$. The value of the factor is less than 0.3 for the step-size $\eta = 1/2$. Therefore the sequence converges to 0, regardless of the initial point.

Note moreover that 0 is always a fixed point of $\bar{\Phi}_\eta$, and that $|\bar{\Phi}_\eta)'(0)| < 1$ for any $\eta \in [7.6, 8.4]$. In that case, 0 is an attracting fixed point, proving the claim. Simulations show that $-0.1$ belongs to the basin of attraction of that point.

$\eta = 2$. For $\eta = 2$, we have $|\bar{\Phi}_\eta(x)| > 2|x|$ for any $x \in \mathbb{R} \setminus \{0\}$ (see Figure 2); therefore the trajectory diverges exponentially fast to $\infty$ from any initial point.

<figure><img src="image_4.png" /><figcaption>Figure 2: Evidence that $\bar{\Phi}_\eta(x) > 2|x|$ for all $x \in \mathbb{R}$, when $\eta = 2$.</figcaption></figure>

$\eta = 3.9$. The map $\psi = (\bar{\Phi}_\eta)^{(4)}$ admits multiple fixed points; see Figure 3(a). For example, Newton's method for computing fixed points (i.e., the sequence defined by $x_{n+1} = x_n - \psi(x_n)/\psi'(x_n)$), initialized at $x = -0.1$, converges numerically to a point $x_\infty \approx -1.35$; its orbit under $\bar{\Phi}_\eta$ is the cycle with approximate values $\{-1.35, 5.92, -1.57, 7.04\}$. The product of $(\bar{\Phi}_\eta)'(u)$ over $u$ in this cycle is approximately $-0.26$, with absolute value less than 1. Therefore, by the next result, this cycle is asymptotically stable for $\bar{\Phi}_\eta$.

**Theorem 5.1** (Prop. 20 Chap. 5 in Block and Coppel [2006]). Let $P = \{x_1, \dots, x_p\}$ be a periodic orbit of a map $\Phi$, and assume that $\Phi$ is differentiable at all $x_i \in P$. If $|\Phi'(x_1) \dots \Phi'(x_p)| < 1$, then $P$ is asymptotically stable.
We illustrate the convergence behavior of this example in Figure 4. Similarly to the previous example, we observe that gradient descent converges to 0 for small step-sizes (in Figure 4(b), for $\eta \le 0.45$), then diverges to $\infty$ for larger step-sizes (at $\eta = 0.5$). An interesting phenomenon occurs when we further increase the step size: the gradient descent iterates start to fill a star-shaped set centered at the origin, (see Figure 4(a)). That is, a set $\mathcal{S}$ such that for any point $x \in \mathcal{S}$, the line segment $[x, 0]$ is contained within the set $\mathcal{S}$. These star-shaped sets appear within a relatively small range of step-size choices; outside this range, we observe divergence.

<figure><img src="image_8.png" /><figcaption>Figure 4: Illustration of Example 4.1. The initial $x_0$ are sampled uniformly at random from the set $[-500, 500]^2$. For Figure 4(a), we sample $n = 1000$ initial $x_0^{(i)}$, $i \in [n]$. The plot shows $x_{n/2}^{(i)}, \dots, x_n^{(i)}$ for each sample $i \in [n]$. Figure 4(b) is based on the $n = 100$ samples for the initial $x_0^{(i)}$, $i \in [n]$. The plot shows the average convergence over these $n$ samples $\frac{1}{n} \sum_{i=1}^n \|x_t^{(i)}\|$.</figcaption></figure>

A complete theoretical analysis of these observations is beyond the scope of this paper. However, we provide a link to existing results in the context of linear iterated function systems (IFS) by Leśniak et al. [2022] and Barnsley and Vince [2011]. To relate our observations and these theoretical results, we need to introduce some notation and terminology. As before, let $\Phi_{\eta,i}: \mathbb{R}^d \to \mathbb{R}^d$ denote the gradient descent operator. Define the IFS as $\mathcal{L} = (\mathbb{R}^d; \Phi_{\eta,i}, i \in [k])$ and its one-parametric version $\mathcal{L}_{\tau} = (\mathbb{R}^d; \frac{1}{\tau}\Phi_{\eta,i}, i \in [k])$. That is, the IFS $\mathcal{L}$ is a collection of $k$ functions $\Phi_{\eta,i}$ from a (metric) set $\mathbb{R}^d$ onto itself. We let $\mathcal{F}^{\Phi} = \{\Phi_{\eta,i}, i \in [k]\}$ denote the set of all functions $\Phi_{\eta,i}$. Let $\mathcal{R}$ denote the set of all non-empty compact subsets of $\mathbb{R}^d$. Recall the Hutchinson operator: $\mathcal{F}: \mathcal{R} \to \mathcal{R}$, with

$$ \mathcal{F}(\mathcal{B}) = \bigcup_{i \in [k]} \Phi_{\eta,i}(\mathcal{B}). $$

A non-empty set $A \subset \mathbb{R}^d$ is an attractor of the IFS if

$$ \mathcal{F}(A) = A \quad \text{and} \quad \forall \mathcal{B} \in \mathcal{R} \quad \lim_{t \to \infty} \mathcal{F}^{(t)}(\mathcal{B}) = A, $$

where $\mathcal{F}^{(t)}$ denotes the $t$-fold composition of the Hutchinson operator $\mathcal{F}$. The limit is taken with respect to the Hausdorff distance. Note that the limit set of the GD operator $\bar{\Phi}_{\eta}$ (recall (9)) for periodic time-varying functions is a subset of the attractor of the IFS $\mathcal{L} = (\mathbb{R}^d; \Phi_{\eta,i}, i \in [k])$. Further, denote by $\rho(\mathcal{F}^{\Phi})$, the joint spectral radius of the IFS (see e.g. Rota and Strang [1960]). In Barnsley and Vince [2011, Theorem 2] and Leśniak et al. [2022], it is shown that an irreducible one-parametric linear IFS $\mathcal{L}_{\tau}$ has a star-shaped, central-symmetric attractor if $\tau = \rho(\mathcal{F}^{\Phi})$. We want to emphasize that these results only apply to linear IFS, while our observations are due to non-linear IFS. Hence, the theory provided in Barnsley and Vince [2011] and Leśniak et al. [2022] does not explain our observation, however, it indicates an interesting future research direction and reveals an interesting link between time-varying optimization problems and IFS.
This shows that there exists an asymptotically stable 4-cycle for $\bar{\Phi}_{\eta}$, i.e., that the sequence $(x_{2n+1})$ will converge to that cycle when initialized in the appropriate (open) set. The result for the whole sequence follows by noting that $x_{2n+2} = \Phi_{\eta,1}(x_{2n+1})$ for any $n \in \mathbb{N}$, and that $\Phi_{\eta,1}$ is continuous.

<figure><img src="image_5.png" /><figcaption>(a) Plot of $I - (\bar{\Phi}_{\eta})^{(4)}$ for $\eta = 3.9$</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>(b) Plot of $I - (\bar{\Phi}_{\eta})^{(3)}$ for $\eta = 6.1$</figcaption></figure>

<caption>Figure 3: Evidence of the existence of periodic points for $\bar{\Phi}_{\eta}$. Each intersection of the plotted curves with the x-axis corresponds to a fixed point of the iterated function, and thus of a periodic orbit for $\bar{\Phi}_{\eta}$. A fixed point of $\Phi_{\eta}^{(3)}$ that is not a fixed point of $\Phi_{\eta}$ has period 3; the existence of such points implies chaos.</caption>

$\eta = 6.1$. The map $\bar{\Phi}_{\eta}$ is an interval map on $[-2.5, 2.5]$. There exists a 3-periodic point $x \approx 0.20$; its orbit is approximately $\{-0.10, 0.20, 0.04\}$. (See Figure 3(b).) Therefore, by the Period-Three theorem [Li and Yorke, 1975], there exists an uncountable scrambled set $S$ in $[-2.5, 2.5]$ for $\bar{\Phi}_{\eta}$. This corresponds to the odd times of the gradient descent trajectory. The result for even times follows again by considering the fact $x_{2n+2} = \Phi_{\eta,1}(x_{2n+1})$, and that $\Phi_{\eta,1}$ is continuous. ■

## 5.4 Iterated function systems and star-shaped attractors.

Another interesting phenomenon we observe for time-varying periodic functions is the convergence to a star-shaped set.

**Example 5.1** (Periodic time-varying problems with star-shaped limit sets). We consider again a 2-periodic sequence as defined in (10), here with dimension $d = 2$ and

$$A_1 = \begin{bmatrix} 3/4 & 0 \\ 0 & 5 \end{bmatrix}, \quad A_2 = \begin{bmatrix} 5 & 1 \\ 1 & 3/4 \end{bmatrix}.$$

We illustrate the convergence behavior of this example in Figure 4. Similarly to the previous example, we observe that gradient descent converges to 0 for small step-sizes (in Figure 4(b), for $\eta \le 0.45$), then diverges to $\infty$ for larger step-sizes (at $\eta = 0.5$). An interesting phenomenon occurs when we further increase the step size: the gradient descent iterates start to fill a star-shaped set centered at the origin, (see Figure 4(a)). That is, a set $S$ such that for any point $x \in S$, the line segment $[x, 0]$ is contained within the set $S$. These star-shaped sets appear within a relatively small range of step-size choices; outside this range, we observe divergence.

A complete theoretical analysis of these observations is beyond the scope of this paper. However, we provide a link to existing results in the context of linear iterated function systems (IFS) by Leśniak et al. [2022] and Barnsley and Vince [2011]. To relate our observations and these theoretical results, we need to introduce some notation and terminology. As before, let $\Phi_{\eta,i}: \mathbb{R}^d \to \mathbb{R}^d$ denote the gradient descent operator. Define the IFS as $\mathcal{L} = (\mathbb{R}^d; \Phi_{\eta,i}, i \in [k])$ and its one-parametric version $\mathcal{L}_{\tau} = (\mathbb{R}^d; \frac{1}{\tau}\Phi_{\eta,i}, i \in [k])$. That is, the IFS $\mathcal{L}$ is a collection of $k$ functions $\Phi_{\eta,i}$ from a (metric) set $\mathbb{R}^d$ onto itself. We let $\mathcal{F}^{\Phi} = \{\Phi_{\eta,i}, i \in [k]\}$ denote the set of all functions $\Phi_{\eta,i}$. Let $\mathcal{R}$ denote the set of all non-empty compact subsets of $\mathbb{R}^d$. Recall the Hutchinson operator: $\mathcal{F}: \mathcal{R} \to \mathcal{R}$, with

$$\mathcal{F}(\mathcal{B}) = \bigcup_{i \in [k]} \Phi_{\eta,i}(\mathcal{B}).$$
## Acknowledgments.

We thank Tim van Erven for his helpful and insightful comments on an earlier version of this manuscript. We also thank Panayotis Mertikopoulos for his feedback on the related literature for Section 2. Further, we thank Mathias Staudigl for some discussions on time-varying games. Sachs, was supported by the Netherlands Organization for Scientific Research (NWO) under grant number VI.Vidi.192.095. Guzmán's research was partially supported by INRIA Associate Teams project, ANID FONDECYT 1210362 grant and ANID Anillo ACT210005 grant.

## References

* *incrementalRegressionLinear*: Linear regression model for incremental learning, 2020. URL https://nl.mathworks.com/help/stats/incrementalregressionlinear.html#refsect-extended-capabilities.

* Ioannis Anagnostides, Ioannis Panageas, Gabriele Farina, and Tuomas Sandholm. On the convergence of no-regret learning dynamics in time-varying games. *Advances in Neural Information Processing Systems*, 36, 2024.

* Dheeraj Baby and Yu-Xiang Wang. Optimal dynamic regret in proper online learning with strongly convex losses and beyond. In *International Conference on Artificial Intelligence and Statistics*, pages 1805–1845. PMLR, 2022.

* Francis Bach and Kfir Y Levy. A universal algorithm for variational inequalities adaptive to smoothness and noise. In *Conference on learning theory*, pages 164–194. PMLR, 2019.

* Michael Barnsley and Andrew Vince. The eigenvalue problem for linear and affine iterated function systems. *Linear Algebra and its Applications*, 435(12):3124–3138, 2011.

* Peter L Bartlett, Elad Hazan, and Alexander Rakhlin. Adaptive online gradient descent. In *Proceedings of the 20th International Conference on Neural Information Processing Systems*, pages 65–72, 2007.

* Heinz H. Bauschke and Patrick L. Combettes. Convex analysis and monotone operator theory in hilbert spaces. *CMS Books in Mathematics*, 2011.

* Heinz H. Bauschke, Sarah M. Moffat, and Xianfu Wang. Firmly nonexpansive mappings and maximally monotone operators: Correspondence and duality. *Set-Valued and Variational Analysis*, 20(1):131–153, July 2011.

* Omar Besbes, Yonatan Gur, and Assaf Zeevi. Non-stationary stochastic optimization. *Operations Research*, 63(5): 1227–1244, 2015.

* Jakub Bielawski, Thiparat Chotibut, Fryderyk Falniowski, Grzegorz Kosiorowski, Michał Misiurewicz, and Georgios Piliouras. Follow-the-regularized-leader routes to chaos in routing games. In *Proceedings of the 38th International Conference on Machine Learning*, volume 139, pages 925–935, 18–24 Jul 2021.

* Louis S Block and William A Coppel. *Dynamics in one dimension*. Springer, 2006.

* Nicolò Cesa-Bianchi and Gábor Lugosi. *Prediction, Learning, and Games*. Cambridge University Press, 2006.

* Xuxing Chen, Krishnakumar Balasubramanian, Promit Ghosal, and Bhavya Agrawalla. From stability to chaos: Analyzing gradient descent dynamics in quadratic regression. *ArXiv*, abs/2310.01687, 2023.

* Steven De Rooij, Tim Van Erven, Peter D. Grünwald, and Wouter M. Koolen. Follow the leader if you can, hedge if you must. *J. Mach. Learn. Res.*, 15(1):1281–1316, jan 2014.

* Benoit Duvocelle, Panayotis Mertikopoulos, Mathias Staudigl, and Dries Vermeulen. Multiagent online learning in time-varying games. *Mathematics of Operations Research*, 48(2):914–941, may 2023.

* Francisco Facchinei and Jong-Shi Pang. *Finite-dimensional variational inequalities and complementarity problems*. Springer, 2003.

* Kenneth Falconer. *Fractal geometry - mathematical foundations and applications*. Wiley, 1990.

* Fryderyk Falniowski and Panayotis Mertikopoulos. On the discrete-time origins of the replicator dynamics: From convergence to instability and chaos. 02 2024.

* Yi Feng, Hu Fu, Qun Hu, Ping Li, Ioannis Panageas, Bo Peng, and Xiao Wang. On the last-iterate convergence in time-varying zero-sum games: Extra gradient succeeds where optimism fails. In *Advances in Neural Information Processing Systems*, volume 36, pages 21933–21944, 2023.

* Tanner Fiez, Ryann Sim, Stratis Skoulakis, Georgios Piliouras, and Lillian Ratliff. Online learning in periodic zero-sum games. In *Advances in Neural Information Processing Systems*, volume 34, pages 10313–10325, 2021.

* Dean P. Foster. Prediction in the worst case. *The Annals of Statistics*, 19(2):1084–1090, April 1991.
<figure><img src="image_7.png" /><figcaption>Figure 4: Illustration of Example 5.1. The initial $x_0$ are sampled uniformly at random from the set $[-500, 500]^2$. For Figure 4(a), we sample $n = 1000$ initial $x_0^{(i)}$, $i \in [n]$. The plot shows $x_{n/2}^{(i)}, \dots, x_n^{(i)}$ for each sample $i \in [n]$. Figure 4(b) is based on the $n = 100$ samples for the initial $x_0^{(i)}$, $i \in [n]$. The plot shows the average convergence over these $n$ samples $\frac{1}{n} \sum_{i=1}^n \|x_t^{(i)}\|$.</figcaption></figure>

A non-empty set $A \subset \mathbb{R}^d$ is an attractor of the IFS if

$$ \mathcal{F}(A) = A \quad \text{and} \quad \forall \mathcal{B} \in \mathcal{R} \quad \lim_{t \to \infty} \mathcal{F}^{(t)}(\mathcal{B}) = A, $$

where $\mathcal{F}^{(t)}$ denotes the $t$-fold composition of the Hutchinson operator $\mathcal{F}$. The limit is taken with respect to the Hausdorff distance. Note that the limit set of the GD operator $\bar{\Phi}_{\eta}$ (recall (9)) for periodic time-varying functions is a subset of the attractor of the IFS $\mathcal{L} = (\mathbb{R}^d; \Phi_{\eta, i}, i \in [k])$. Further, denote by $\rho(\mathcal{F}^{\Phi})$, the joint spectral radius of the IFS (see e.g. Rota and Strang [1960]). In Barnsley and Vince [2011, Theorem 2] and Leśniak et al. [2022], it is shown that an irreducible one-parametric linear IFS $\mathcal{L}_{\tau}$ has a star-shaped, central-symmetric attractor if $\tau = \rho(\mathcal{F}^{\Phi})$. We want to emphasize that these results only apply to linear IFS, while our observations are due to non-linear IFS. Hence, the theory provided in Barnsley and Vince [2011] and Leśniak et al. [2022] does not explain our observation, however, it indicates an interesting future research direction and reveals an interesting link between time-varying optimization problems and IFS.

# 6 Discussion, Limitations and Future Work

We studied the problem of tracking equilibria in time-varying variational inequality problems. Positive results include sublinear tracking error bounds for tame instances and contractive algorithms. These bounds are proven to be tight for this class of algorithms. Other positive results concern the case of periodic instances with strongly monotone operators, where we propose meta-algorithms for aggregating multiple experts. These algorithms result on logarithmic or even constant tracking bounds. Finally, we explore gradient-type methods with constant stepsize for periodic instances, as an agnostic approach to tackle periodic instances with unknown periods. Our results show an intricate landscape where convergence, chaos, and nonconvergence may arise, making the use of these algorithms difficult even in simple cases.

We discuss next some of the limitations in this work. The algorithm of Section 3.1 is computationally more efficient than the algorithm in Section 3.2. However, in contrast to the latter, it requires a bounded domain and provides (only) logarithmic guarantees. It would be an interesting question for future work to combine the merits of both approaches. Another limitation becomes clear when considering the special case of games: In this case, the aggregation schemes discussed in this section may be unnatural, as they require coordination between the players. Decoupling the aggregation in general cases would be an interesting direction for further investigation.

We also remark that both results in Section 3.1 and 3.2 require deterministic feedback. An interesting direction is to consider more robust models, in which we would observe the inexact values of the operators. Examples of such models could be a stochastic feedback model or time-varying VI problems that admit periodic behavior only in the
E. Hazan. *Introduction to Online Convex Optimization, second edition*. Adaptive Computation and Machine Learning series. MIT Press, 2022. ISBN 9780262370127.

Liam Hodgkinson, Umut Simsekli, Rajiv Khanna, and Michael W. Mahoney. Generalization properties of stochastic optimizers via trajectory analysis. *ArXiv*, abs/2108.00781, 2021.

Anatoli B. Juditsky and Arkadi S. Nemirovsky. Signal recovery by stochastic optimization. *Autom. Remote. Control.*, 80(10):1878–1893, 2019.

Hamed Karimi, Julie Nutini, and Mark W. Schmidt. Linear convergence of gradient and proximal-gradient methods under the polyak-lojasiewicz condition. In *ECML/PKDD*, 2016.

F P Kelly, A K Maulloo, and D K H Tan. Rate control for communication networks: shadow prices, proportional fairness and stability. *Journal of the Operational Research Society*, 49(3):237–252, March 1998.

David Kinderlehrer and Guido Stampacchia. *An Introduction to Variational Inequalities and Their Applications*. SIAM, January 2000. ISBN 9780898719451.

R Klén, M Visuri, and M Vuorinen. On Jordan type inequalities for hyperbolic functions. *Journal of Inequalities and Applications*, 2010(1):362548, 2010.

Galina M Korpelevich. The extrag gradient method for finding saddle points and other problems. *Matecon*, 12:747–756, 1976.

Krzysztof Leśniak, Nina Snigireva, Filip Strobin, and Andrew Vince. Transition phenomena for the attractor of an iterated function system*. *Nonlinearity*, 35(10):5396–5426, September 2022.

Tien-Yien Li and James A Yorke. Period three implies chaos. *The American Mathematical Monthly*, 82(10):985–992, 1975.

Konstantin Mishchenko, Ahmed Khaled, and Peter Richtarik. Random reshuffling: Simple analysis with vast improvements. In *Advances in Neural Information Processing Systems*, volume 33, pages 17309–17320, 2020.

Aryan Mokhtari, Shahin Shahrampour, Ali Jadbabaie, and Alejandro Ribeiro. Online optimization in dynamic environments: Improved regret rates for strongly convex problems. In *2016 IEEE 55th Conference on Decision and Control (CDC)*, pages 7195–7201. IEEE, 2016.

Angelia Nedic and Dimitri P. Bertsekas. Incremental subgradient methods for nondifferentiable optimization. *SIAM Journal on Optimization*, 12(1):109–138, 2001.

L.M. Nguyen, Q. Tran-Dinh, D.T. Phan, P.H. Nguyen, and M. van Dijk. A unified convergence analysis for shuffling-type gradient methods. *Journal of Machine Learning Research*, 22, 2021. ISSN 1532-4435.

Georgios Piliouras and Fang-Yi Yu. Multi-agent performative prediction: From global stability and optimality to chaos. In *Proceedings of the 24th ACM Conference on Economics and Computation*, EC '23. ACM, July 2023.

Alexander Rakhlin and Karthik Sridharan. Online learning with predictable sequences. In *Conference on Learning Theory*, pages 993–1019. PMLR, 2013.

Adrian Rivera Cardoso, Jacob Abernethy, He Wang, and Huan Xu. Competing against Nash equilibria in adversarially changing zero-sum games. In *Proceedings of the 36th International Conference on Machine Learning*, volume 97 of *Proceedings of Machine Learning Research*, pages 921–930, 2019.

R. Tyrrell Rockafellar and Roger J. B. Wets. *Variational Analysis*. Springer, 1998.

J. B. Rosen. Existence and uniqueness of equilibrium points for concave n-person games. *Econometrica*, 33(3): 520–534, 1965. ISSN 00129682, 14680262.

Gian-Carlo Rota and W. Gilbert Strang. A note on the joint spectral radius. 1960.

Sylvie Ruette. *Chaos on the Interval*, volume 67. American Mathematical Soc., 2017.

Shai Shalev-Shwartz and Shai Ben-David. *Understanding Machine Learning: From Theory to Algorithms*. Cambridge University Press, USA, 2014. ISBN 1107057132.

Tim Van Erven, Wouter M Koolen, and Dirk Van der Hoeven. Metagrad: Adaptation using multiple learning rates in online learning. *Journal of Machine Learning Research*, 22(161):1–61, 2021.

Vladimir G Vovk. A game of prediction with expert advice. In *Proceedings of the eighth annual conference on Computational learning theory*, pages 51–60, 1995.

Chen-Yu Wei, Chung-Wei Lee, Mengxiao Zhang, and Haipeng Luo. Linear last-iterate convergence in constrained saddle-point optimization. published at ICLR'21, 2021.

Yu-Hu Yan, Peng Zhao, and Zhi-Hua Zhou. Fast rates in time-varying strongly monotone games. In *Proceedings of the 40th International Conference on Machine Learning*, volume 202, pages 39138–39164, 23–29 Jul 2023.
limit. These feedback models would capture application examples, such as fluctuating markets, in a more realistic setting. Inexact feedback models are inherently more challenging; the convergence guarantees in this setting are weaker, and the aggregation step would be more difficult to implement.

Another interesting future research direction emerges from the relation of time-varying games to stochastic games. Stochastic games are a broad class of time-varying games, where the state of the game evolves depending on both the players' actions as well as random events. This setting contains both tame time-varying games and periodic ones, but it is a much broader class. Understanding the class of stochastic games for which tracking is attainable is then a very interesting and important direction for future research.

## Acknowledgments.

We thank Tim van Erven for his helpful and insightful comments on an earlier version of this manuscript. We also thank Panayotis Mertikopoulos for his feedback on the related literature for Section 2. Further, we thank Mathias Staudigl for some discussions on time-varying games. Work was done while Sachs was supported by the Netherlands Organization for Scientific Research (NWO) under grant number VI.Vidi.192.095, now Sachs is supported by the European Union - Next Generation EU, PRIN 2022, CUP:J53D23007170001. Guzmán's research was partially supported by INRIA Associate Teams project, ANID FONDECYT 1210362 grant and ANID Anillo ACT210005 grant.

## References

* *incrementalRegressionLinear*: Linear regression model for incremental learning, 2020. URL https://nl.mathworks.com/help/stats/incrementalregressionlinear.html#refsect-extended-capabilities.

* I. Anagnostides, I. Panageas, G. Farina, and T. Sandholm. On the convergence of no-regret learning dynamics in time-varying games. *Advances in Neural Information Processing Systems*, 36, 2024.

* K. Antonakopoulos, D. Q. Vu, V. Cevher, K. Levy, and P. Mertikopoulos. Undergrad: A universal black-box optimization method with almost dimension-free convergence rate guarantees. In *International Conference on Machine Learning*, pages 772–795. PMLR, 2022.

* L. Armijo. Minimization of functions having lipschitz continuous first partial derivatives. *Pacific Journal of mathematics*, 16(1):1–3, 1966.

* D. Baby and Y.-X. Wang. Optimal dynamic regret in proper online learning with strongly convex losses and beyond. In *International Conference on Artificial Intelligence and Statistics*, pages 1805–1845. PMLR, 2022.

* M. Barnsley and A. Vince. The eigenvalue problem for linear and affine iterated function systems. *Linear Algebra and its Applications*, 435(12):3124–3138, 2011.

* P. L. Bartlett, E. Hazan, and A. Rakhlin. Adaptive online gradient descent. In *Proceedings of the 20th International Conference on Neural Information Processing Systems*, pages 65–72, 2007.

* L. Baudin, M. Scarsini, and X. Venel. Strategic behavior and no-regret learning in queueing systems. *ArXiv*, abs/2302.03614, 2023.

* H. H. Bauschke and P. L. Combettes. Convex analysis and monotone operator theory in hilbert spaces. *CMS Books in Mathematics*, 2011.

* H. H. Bauschke, S. M. Moffat, and X. Wang. Firmly nonexpansive mappings and maximally monotone operators: Correspondence and duality. *Set-Valued and Variational Analysis*, 20(1):131–153, July 2011.

* O. Besbes, Y. Gur, and A. Zeevi. Non-stationary stochastic optimization. *Operations Research*, 63(5):1227–1244, 2015.

* J. Bielawski, T. Chotibut, F. Falniowski, G. Kosiorowski, M. Misiurewicz, and G. Piliouras. Follow-the-regularized-leader routes to chaos in routing games. In *Proceedings of the 38th International Conference on Machine Learning*, volume 139, pages 925–935, 18–24 Jul 2021.

* L. S. Block and W. A. Coppel. *Dynamics in one dimension*. Springer, 2006.

* N. Cesa-Bianchi and G. Lugosi. *Prediction, Learning, and Games*. Cambridge University Press, 2006.

* X. Chen, K. Balasubramanian, P. Ghosal, and B. Agrawalla. From stability to chaos: Analyzing gradient descent dynamics in quadratic regression. *ArXiv*, abs/2310.01687, 2023.

* J. M. Cohen, S. Kaur, Y. Li, J. Z. Kolter, and A. Talwalkar. Gradient descent on neural networks typically occurs at the edge of stability. *ArXiv*, abs/2103.00065, 2021.
Lijun Zhang, Shiyin Lu, and Zhi-Hua Zhou. Adaptive online learning in dynamic environments. In Proceedings of the 32nd International Conference on Neural Information Processing Systems, NIPS'18, pages 1330–1340, 2018.

Mengxiao Zhang, Peng Zhao, Haipeng Luo, and Zhi-Hua Zhou. No-regret learning in time-varying zero-sum games. In International Conference on Machine Learning, pages 26772–26808. PMLR, 2022.

Peng Zhao and Lijun Zhang. Improved analysis for dynamic regret of strongly convex and smooth functions. In Learning for Dynamics and Control, pages 48–59. PMLR, 2021.

Martin Zinkevich. Online convex programming and generalized infinitesimal gradient ascent. In Proceedings of the 20th international conference on machine learning (icml-03), pages 928–936, 2003.

# Appendix

## A Proofs for Section 1.

**Lemma A.1 (Lower Bound).** For any algorithm A, there exists a sequence $(F_t)_{t \in [T]}$ such that $\tau(A) \ge \frac{1}{4} \sum_{t=1}^T \|Z_t^* - Z_{t+1}^*\|^2$ and $\sum_{t=1}^T \|Z_t^* - Z_{t+1}^*\|^2 = \frac{1}{4} D^2 T$.

*Proof.* We show this lower bound by defining a sequence $(Z_t^*)_{t \in [T]}$ with $Z_t^* \in \mathcal{Z}$. The claim follows by choosing $F_t$ as the gradient of a convex function $f_t$ with a unique minimizer at $Z_t^*$. This implies that $Z_t^*$ solves the (VIP). We assume that $\mathcal{Z} = [-1, 1] \subseteq \mathbb{R}$ and choose $Z_t^* \in \{-1, 0, 1\}$. For any $t \in \mathbb{N}$, assume wlog $Z_t \ge 0$. We distinguish two cases: (1) $Z_{t-1}^* \in \{0, 1\}$ and (2) $Z_{t-1}^* = -1$. For case (1), we set $Z_t^* = -1$ and directly obtain $\|Z_t - Z_t^*\|^2 \ge \frac{1}{4} \|Z_{t-1}^* - Z_t^*\|^2$. For case (2), if $Z_t \ge \frac{1}{2}$, we set $Z_t^* = 0$ otherwise $Z_t^* = 1$. In both cases, we directly obtain $\|Z_t - Z_t^*\|^2 \ge \frac{1}{4} \|Z_{t-1}^* - Z_t^*\|^2$. ■

## B Proofs for Section 2.

### B.1 Proofs: Example 2.2

We adopt an example by Karimi et al. [2016]. Recall the time-varying game from (2) with $a \in [0, 1]$

$$ \ell(x, y) = x^2 + 3 \sin^2(x) + a \sin^2(x) \sin^2(y) - y^2 - 3 \sin^2(y). $$

This game is not convex-concave. Its gradients with respect to $x$ are for any $x, y$,

$$ \nabla_x \ell(x, y) = x + 3 \sin(x) \cos(x) + a \sin(x) \cos(x) \sin^2(y) = 2x + \sin(2x)(3 + a \sin^2(y))/2, $$

$$ \nabla_y \ell(x, y) = -y - 3 \sin(y) \cos(y) + a \sin(y) \cos(y) \sin^2(x) = -2y - \sin(2x)(3 - a \sin^2(x))/2. $$

We now show that the game defined in (2) is

1. not a monotone game,

2. satisfies the Restricted Secant Inequality (RSI) in $x$ for any fixed $y$ (and in $y$ for any fixed $x$, respectively), i.e., there exists a $\mu_{\text{RSI}} > 0$ such that

$$ \forall x \in \mathbb{R}: \langle \nabla_x \ell(x, y), x - x^* \rangle \ge \mu_{\text{RSI}} \|x - x^*\|^2, $$

where $x^*$ denotes the projection of $x$ onto the set of minimizers.

3. GDA satisfies the contraction assumption.

**Proposition B.1.** For any $x, y \in \mathbb{R}$, the functions $\ell(\cdot, y)$ and $\ell(x, \cdot)$ satisfy the RSI with $\mu = 1/4$.

*Proof.* We want to show that

$$ (2x + \sin(2x)(3 + a \sin^2(y))/2)(x - x^*) \ge \mu(x - x^*)^2. \quad (11) $$

First note that for any $y$, $\ell(x, y)$ has a unique minimizer $x^* = 0$. Further, for $x = 0$, (11) holds for any $\mu$. Hence, assume $x \ne 0$. Then (11) is satisfied if

$$ 2 + \frac{\sin(2x)}{2x}(3 + a \sin^2(y)) \ge \mu, $$
S. De Rooij, T. Van Erven, P. D. Grünwald, and W. M. Koolen. Follow the leader if you can, hedge if you must. *J. Mach. Learn. Res.*, 15(1):1281–1316, jan 2014.

Y. Deng, V. Mirrokni, and S. Zuo. Non-clairvoyant dynamic mechanism design with budget constraints and beyond. In *Proceedings of the 22nd ACM Conference on Economics and Computation*, pages 369–369, 2021.

B. Duvocelle, P. Mertikopoulos, M. Staudigl, and D. Vermeulen. Multiagent online learning in time-varying games. *Mathematics of Operations Research*, 48(2):914–941, may 2023.

F. Facchinei and C. Kanzow. Generalized Nash equilibrium problems. *4OR*, 5(3):173–210, sep 2007.

F. Facchinei and J.-S. Pang. *Finite-dimensional variational inequalities and complementarity problems*. Springer, 2003.

K. Falconer. *Fractal geometry - mathematical foundations and applications*. Wiley, 1990.

F. Falniowski and P. Mertikopoulos. On the discrete-time origins of the replicator dynamics: From convergence to instability and chaos. 02 2024.

Y. Feng, H. Fu, Q. Hu, P. Li, I. Panageas, B. Peng, and X. Wang. On the last-iterate convergence in time-varying zero-sum games: Extra gradient succeeds where optimism fails. In *Advances in Neural Information Processing Systems*, volume 36, pages 21933–21944, 2023.

T. Fiez, R. Sim, S. Skoulakis, G. Piliouras, and L. Ratliff. Online learning in periodic zero-sum games. In *Advances in Neural Information Processing Systems*, volume 34, pages 10313–10325, 2021.

D. P. Foster. Prediction in the worst case. *The Annals of Statistics*, 19(2):1084–1090, April 1991.

E. Hazan. *Introduction to Online Convex Optimization, second edition*. Adaptive Computation and Machine Learning series. MIT Press, 2022. ISBN 9780262370127.

L. Hodgkinson, U. Simsekli, R. Khanna, and M. W. Mahoney. Generalization properties of stochastic optimizers via trajectory analysis. *ArXiv*, abs/2108.00781, 2021.

A. B. Juditsky and A. S. Nemirovsky. Signal recovery by stochastic optimization. *Autom. Remote. Control.*, 80(10): 1878–1893, 2019.

H. Karimi, J. Nutini, and M. W. Schmidt. Linear convergence of gradient and proximal-gradient methods under the polyak-lojasiewicz condition. In *ECML/PKDD*, 2016.

A. Kavis, K. Y. Levy, F. Bach, and V. Cevher. Unixgrad: A universal, adaptive algorithm with optimal guarantees for constrained optimization. *Advances in neural information processing systems*, 32, 2019.

F. P. Kelly, A. K. Maulloo, and D. K. H. Tan. Rate control for communication networks: shadow prices, proportional fairness and stability. *Journal of the Operational Research Society*, 49(3):237–252, Mar. 1998.

D. Kinderlehrer and G. Stampacchia. *An Introduction to Variational Inequalities and Their Applications*. SIAM, Jan. 2000. ISBN 9780898719451.

R. Klén, M. Visuri, and M. Vuorinen. On Jordan type inequalities for hyperbolic functions. *Journal of Inequalities and Applications*, 2010(1):362548, 2010.

G. M. Korpelevich. The extragrwardient method for finding saddle points and other problems. *Matecon*, 12:747–756, 1976.

T. Larsson and M. Patriksson. A class of gap functions for variational inequalities. *Mathematical Programming*, 64: 53–79, 1994.

K. Leśniak, N. Snigireva, F. Strobin, and A. Vince. Transition phenomena for the attractor of an iterated function system*. *Nonlinearity*, 35(10):5396–5426, Sept. 2022.

K. Y. Levy, A. Yurtsever, and V. Cevher. Online adaptive methods, universality and acceleration. *Advances in neural information processing systems*, 31, 2018.

T.-Y. Li and J. A. Yorke. Period three implies chaos. *The American Mathematical Monthly*, 82(10):985–992, 1975.

P. Mertikopoulos and M. Staudigl. Equilibrium tracking and convergence in dynamic games. In *2021 60th IEEE Conference on Decision and Control (CDC)*, pages 930–935, United States, 2021. IEEE.

K. Mishchenko, A. Khaled, and P. Richtarik. Random reshuffling: Simple analysis with vast improvements. In *Advances in Neural Information Processing Systems*, volume 33, pages 17309–17320, 2020.

A. Mokhtari, S. Shahrampour, A. Jadbabaie, and A. Ribeiro. Online optimization in dynamic environments: Improved regret rates for strongly convex problems. In *2016 IEEE 55th Conference on Decision and Control (CDC)*, pages 7195–7201. IEEE, 2016.
Using the fact that $\sin(x) \in [-1, 1]$ and that the minimum of $\sin(x)/x$ is attained at $\{-\frac{3}{2}\pi, \frac{3}{2}\pi\}$ (see, e.g., Klén et al. [2010]), $\sin(2x)/(2x) \ge -0.22$. Therefore, since $3 \le (3 + a \sin^2(y)) \le 4$ for any $a$ and $y$, we have

$$ \frac{\sin(2x)}{2x}(3 + a \sin^2(y)) \ge -0.88, \quad \text{so} \quad 2 + \frac{\sin(2x)}{2x}(3 + a \sin^2(y)) \ge 1.12, $$

proving the claim.

<figure><img src="image_9.png" /></figure>

**Proposition B.2.** The game defined in (2) is not convex and a fortiori not strongly convex.

*Proof.* The second derivative with respect to $x$ evaluated at $x, y$ is $2 + \cos(2x)(3 + a \sin^2(y))$. which can be negative, e.g. at $x = \pi/2$ and $y = 0$.

<figure><img src="image_10.png" /></figure>

**Proposition B.3.** Contraction is satisfied for problem (2) with GDA with $C = (1 - \mu_{\text{RSI}}^2/L^2)^{1/2}$ where $\mu_{\text{RSI}}$ denotes the RSI parameter.

*Proof.* As for the strongly-convex case:

$$ \|x_2 - x^*\|^2 = \|x - x^* - \eta \nabla_x \ell(x, y)\|^2 = \|x - x^*\|^2 - 2\eta \langle x - x^*, \nabla_x \ell(x, y) \rangle + \eta^2 \|\nabla_x \ell(x, y)\|^2. $$

By RSI

$$ \langle x - x^*, \nabla_x \ell(x, y) \rangle + \langle y - y^*, -\nabla_y \ell(x, y) \rangle \ge \mu_{\text{RSI}} \|x - x^*\|^2 + \mu_{\text{RSI}} \|y - y^*\|^2, $$

thus,

$$ \|x_2 - x^*\|^2 + \|y_2 - y^*\|^2 \le (1 - 2\eta\mu_{\text{RSI}})(\|x_1 - x^*\|^2 + \|y_1 - y^*\|^2) + \eta^2 (\|\nabla_x \ell(x, y)\|^2 + \|\nabla_y \ell(x, y)\|^2). $$

The result follows as before via L-smoothness.

<figure><img src="image_11.png" /></figure>

## B.2 Tightness.

### Proof of Theorem 2.2.

*Proof.* We now provide an example of a sequence of problems, with a C-contractive algorithm, such that the bound given in Theorem 2.1 is tight. Our example is an elementary analysis of gradient descent on a sequence of 1-dimensional quadratic optimization problems, where the optima drift linearly. Consider a sequence of quadratic problems over $\mathbb{R}$ defined by

$$ f_t(x) = \frac{(x - c_t)^2}{2}; $$

for our example, let us assume that $(c_t)$ forms an arithmetic progression with common difference $-b$, i.e., $c_{t+1} = c_t - b$ for all $t \ge 1$. Consider also the gradient descent algorithm started at $x_1 \in \mathbb{R}$, and assume furthermore that $x_1 - c_1 \ge b/(1-C)$.

Then for any $t \in \mathbb{R}$, gradient descent started at $x_1 \in \mathbb{R}$ with step-size $(1-C) \in (0, 1]$ is a C-contraction since for any $t \ge 1$, the update rule is

$$ x_{t+1} - c_t = x_t - (1-C)(x_t - c_t) - c_t = C(x_t - c_t). $$

Let us introduce the short-hand notation $d_t = x_t - c_t$. Then for any $t \ge 1$,

$$ d_{t+1} = x_{t+1} - c_{t+1} = x_{t+1} - c_t + c_t - c_{t+1} = Cd_t + b $$

The sequence $(d_t)$ therefore forms an arithmetico-geometric progression and so for any $t \ge 1$, using our assumption that $x_1 - c_1 \ge b/(1-C)$,

$$ d_t - \frac{b}{1-C} = C^{t-1} \left( d_1 - \frac{b}{1-C} \right) \ge 0. $$

Therefore, for any $t \ge 1$, the squared-error $d_t^2$ can be lower bounded as

$$ \begin{aligned} d_t^2 &= \left( d_t - \frac{b}{1-C} + \frac{b}{1-C} \right)^2 \\ &= \left( d_t - \frac{b}{1-C} \right)^2 + \frac{b^2}{(1-C)^2} + 2 \left( d_t - \frac{b}{1-C} \right) \frac{b}{1-C} \ge \frac{b^2}{(1-C)^2}. \end{aligned} $$

Conclude by summing the inequality above over $t \in [T]$,

$$ \sum_{t=1}^{T} (z_t - c_t)^2 \ge \frac{Tb^2}{(1-C)^2} \ge \frac{1}{(1-C)^2} \sum_{t=2}^{T} |c_t - c_{t-1}|^2, $$

proving the claim.

<figure><img src="image_12.png" /></figure>
T. Nedelec, C. Calauzènes, N. El Karoui, V. Perchet, et al. Learning in repeated auctions. *Foundations and Trends in Machine Learning*, 15(3):176–334, 2022.

A. Nedic and D. P. Bertsekas. Incremental subgradient methods for nondifferentiable optimization. *SIAM Journal on Optimization*, 12(1):109–138, 2001.

G. Piliouras and F.-Y. Yu. Multi-agent performative prediction: From global stability and optimality to chaos. In *Proceedings of the 24th ACM Conference on Economics and Computation*, EC '23. ACM, July 2023.

J. Rahme, S. Jelassi, and S. M. Weinberg. Auction learning as a two-player game. In *International Conference on Learning Representations*, 2021.

A. Rakhlin and K. Sridharan. Online learning with predictable sequences. In *Conference on Learning Theory*, pages 993–1019. PMLR, 2013.

A. Rivera Cardoso, J. Abernethy, H. Wang, and H. Xu. Competing against Nash equilibria in adversarially changing zero-sum games. In *Proceedings of the 36th International Conference on Machine Learning*, volume 97 of *Proceedings of Machine Learning Research*, pages 921–930, 2019.

R. T. Rockafellar and R. J. B. Wets. *Variational Analysis*. Springer, 1998.

J. B. Rosen. Existence and uniqueness of equilibrium points for concave n-person games. *Econometrica*, 33(3): 520–534, 1965. ISSN 00129682, 14680262.

G.-C. Rota and W. G. Strang. A note on the joint spectral radius. 1960.

S. Ruette. *Chaos on the Interval*, volume 67. American Mathematical Soc., 2017.

S. Shalev-Shwartz and S. Ben-David. *Understanding Machine Learning: From Theory to Algorithms*. Cambridge University Press, USA, 2014. ISBN 1107057132.

T. Van Erven, W. M. Koolen, and D. Van der Hoeven. Metagrad: Adaptation using multiple learning rates in online learning. *Journal of Machine Learning Research*, 22(161):1–61, 2021.

V. G. Vovk. A game of prediction with expert advice. In *Proceedings of the eighth annual conference on Computational learning theory*, pages 51–60, 1995.

C.-Y. Wei, C.-W. Lee, M. Zhang, and H. Luo. Linear last-iterate convergence in constrained saddle-point optimization. published at ICLR'21, 2021.

Y.-H. Yan, P. Zhao, and Z.-H. Zhou. Fast rates in time-varying strongly monotone games. In *Proceedings of the 40th International Conference on Machine Learning*, volume 202, pages 39138–39164, 23–29 Jul 2023.

L. Zhang, S. Lu, and Z.-H. Zhou. Adaptive online learning in dynamic environments. In *Proceedings of the 32nd International Conference on Neural Information Processing Systems*, NIPS'18, pages 1330–1340, 2018.

M. Zhang, P. Zhao, H. Luo, and Z.-H. Zhou. No-regret learning in time-varying zero-sum games. In *International Conference on Machine Learning*, pages 26772–26808. PMLR, 2022.

P. Zhao and L. Zhang. Improved analysis for dynamic regret of strongly convex and smooth functions. In *Learning for Dynamics and Control*, pages 48–59. PMLR, 2021.

M. Zinkevich. Online convex programming and generalized infinitesimal gradient ascent. In *Proceedings of the 20th international conference on machine learning (icml-03)*, pages 928–936, 2003.

A Proofs for Section 1.

**Lemma A.1 (Lower Bound).** For any algorithm $\mathcal{A}$, there exists a sequence $(F_t)_{t \in [T]}$ such that $\tau(\mathcal{A}) \ge \frac{1}{4} \sum_{t=1}^T \|Z_t^* - Z_{t+1}^*\|^2$ and $\sum_{t=1}^T \|Z_t^* - Z_{t+1}^*\|^2 = \frac{1}{4} D^2 T$.

*Proof.* Proof. We show this lower bound by defining a sequence $(Z_t^*)_{t \in [T]}$ with $Z_t^* \in \mathcal{Z}$. The claim follows by choosing $F_t$ as the gradient of a convex function $f_t$ with a unique minimizer at $Z_t^*$. This implies that $Z_t^*$ solves the (VIP). We assume that $\mathcal{Z} = [-1, 1] \subseteq \mathbb{R}$ and choose $Z_t^* \in \{-1, 0, 1\}$. For any $t \in \mathbb{N}$, we distinguish two cases: $Z_t \ge 0$ and $Z_t < 0$. First assume $Z_t \ge 0$. We distinguish again two cases: (1) $Z_{t-1}^* \in \{0, 1\}$ and (2) $Z_{t-1}^* = -1$. For case (1), we set $Z_t^* = -1$ and directly obtain $\|Z_t - Z_t^*\|^2 \ge \frac{1}{4} \|Z_{t-1}^* - Z_t^*\|^2$. For case (2), if $Z_t \ge \frac{1}{2}$, we set $Z_t^* = 0$ otherwise $Z_t^* = 1$. In both cases, we directly obtain $\|Z_t - Z_t^*\|^2 \ge \frac{1}{4} \|Z_{t-1}^* - Z_t^*\|^2$. We build $Z_{t-1}^*$ using a symmetric construction when $Z_t < 0$, which completes the proof. ■
# C Proofs of main positive results.

## C.1 Proof of Theorem 3.1

*Proof.* The setting for aggregation we consider here is an instance of Prediction with Expert Advice, see Cesa-Bianchi and Lugosi [2006, Chapter 2], where the $K$ experts provide the advice $z_t^{(i)}$, the outcomes are the pairs $(F_t(z_t), z_t)$, and where the loss is the function over $\mathcal{Z} \times ((\mathcal{Z})^* \times \mathcal{Z})$:

$$ \mathcal{L}(a; g, z) = \langle g, a \rangle + \frac{\mu}{2} \|a - z\|^2 = \frac{\mu}{2} \|a - z - \frac{1}{\mu}g\|^2 - \frac{1}{\mu} \|g\|^2 + \langle z, g \rangle; $$

we apply it to $(g, z) = (F_t(z_t), z_t)$, thus the predictions will be bounded by $D$ and the outcome vectors by $D + G/\mu$. This loss function, which is an instance of the square loss, is $\lambda$-mixable for the mean (or, in other words, $\lambda$-exp-concave) for $\lambda \le (4\mu(D + G/\mu)^2)^{-1}$, as long as $\|x + (1/\mu)g\| \le D + G/\mu$ and $\|a\| \le D$. Playing the Exponential Weights updates from (6) in this context is equivalent to using Vovk's aggregating algorithm Vovk [1995], which guarantees for $\lambda$-mixable losses that

$$ \sum_{t=1}^{T} \mathcal{L}(Z_t; F_t(Z_t), Z_t) \le \min_{i \in [K]} \sum_{t=1}^{T} \mathcal{L}(Z_t^{(i)}; F_t(Z_t), Z_t) + \frac{\log K}{\lambda}. $$

This can be written as

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t \rangle \le \min_{i \in [K]} \sum_{t=1}^{T} \langle F_t(Z_t), Z_t^{(i)} \rangle + \frac{\mu}{2} \|Z_t^{(i)} - Z_t\|^2 + \frac{\log K}{\lambda}. \quad (12) $$

Now note that for any $u \in \mathcal{Z}$,

$$ \begin{aligned} \langle \tilde{F}_t(Z_t^{(i)}), Z_t^{(i)} - u \rangle - \frac{\mu}{2} \|Z_t^{(i)} - u\|^2 &= \langle F_t(Z_t) + \mu(Z_t^{(i)} - Z_t), Z_t^{(i)} - u \rangle - \frac{\mu}{2} \|Z_t^{(i)} - u\|^2 \\ &= \langle F_t(Z_t), Z_t^{(i)} - u \rangle + \mu \langle Z_t^{(i)} - Z_t, Z_t^{(i)} - u \rangle - \frac{\mu}{2} \|Z_t^{(i)} - u\|^2 \\ &= \langle F_t(Z_t), Z_t^{(i)} - u \rangle - \frac{\mu}{2} \|Z_t - u\|^2 + \frac{\mu}{2} \|Z_t^{(i)} - Z_t\|^2. \end{aligned} $$

Therefore, plugging this bound in (12), we see that for any sequence of comparators $(C_t) \in \mathcal{Z}^T$,

$$ \begin{aligned} \sum_{t=1}^{T} \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 &\le \min_{i \in [K]} \sum_{t=1}^{T} \langle F_t(Z_t), Z_t^{(i)} - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 + \frac{\mu}{2} \|Z_t^{(i)} - Z_t\|^2 + \frac{\log K}{\lambda} \\ &= \min_{i \in [K]} \sum_{t=1}^{T} \langle \tilde{F}_t(Z_t^{(i)}), Z_t^{(i)} - C_t \rangle - \frac{\mu}{2} \|Z_t^{(i)} - C_t\|^2 + \frac{\log K}{\lambda}. \end{aligned} $$

If in addition, $(C_t)$ is $k$-periodic, then, by applying the guarantees from Corollary 3.1 to the operators $\tilde{F}_t$ which are upper bounded in norm by $G + \mu D$, we obtain

$$ \min_{i \in [k]} \sum_{t=1}^{T} \langle \tilde{F}_t(Z_t^{(i)}), Z_t^{(i)} - C_t \rangle - \frac{\mu}{2} \|Z_t^{(i)} - C_t\|^2 \le \frac{k(G + \mu D)^2}{2\mu} \log \left( \frac{T}{k} + 1 \right). $$

Replace $\lambda$ by its value to conclude the proof.

<figure><img src="image_13.png" /></figure>

## C.2 Proof of Theorem 3.2

We separate the proof into three distinct lemmas. We start with an analysis of the aggregation method, showing that the regret of the full procedure can be controlled by the cumulative distance of the correctly tuned base algorithm to the comparators, plus the cumulative mixability gap. We then proceed to bound the said cumulative mixability gap, by showing that the instantaneous mixability gaps are eventually 0, making the sum bounded independently of $T$. (More accurately, the dependence is only in $\hat{D}_T$ and $\hat{G}_T$, which will then be bounded independently.)

Recall that $(x)_+$ denotes the positive part of a real number $x$, and that

$$ \hat{D}_T = \max_{t \in [T], i, j \in [K]} \|Z_t^{(i)} - Z_t^{(j)}\| \quad \text{and} \quad \hat{G}_T = \max_{t \in [T], i, j \in [K]} \|F_t(Z_t^{(i)})\|. $$
# B Proofs for Section 2.

## B.1 Proofs: Example 2.2

We adopt an example by Karimi et al. [2016]. Recall the time-varying game from (1) with $a \in [0, 1]$

$$ \ell(x, y) = x^2 + 3 \sin^2(x) + a \sin^2(x) \sin^2(y) - y^2 - 3 \sin^2(y). $$

This game is not convex-concave. Its gradients with respect to $x$ are for any $x, y$,

$$ \nabla_x \ell(x, y) = x + 3 \sin(x) \cos(x) + a \sin(x) \cos(x) \sin^2(y) = 2x + \sin(2x)(3 + a \sin^2(y))/2, $$

$$ \nabla_y \ell(x, y) = -y - 3 \sin(y) \cos(y) + a \sin(y) \cos(y) \sin^2(x) = -2y - \sin(2x)(3 - a \sin^2(x))/2. $$

We now show that the game defined in (1) is

1. not a monotone game,

2. satisfies the Restricted Secant Inequality (RSI) in $x$ for any fixed $y$ (and in $y$ for any fixed $x$, respectively), i.e., there exists a $\mu_{RSI} > 0$ such that

$$ \forall x \in \mathbb{R}: \langle \nabla_x \ell(x, y), x - x^* \rangle \geq \mu_{\text{RSI}} \|x - x^*\|^2, $$

where $x^*$ denotes the projection of $x$ onto the set of minimizers.

3. GDA satisfies the contraction assumption.

**Proposition B.1.** For any $x, y \in \mathbb{R}$, the functions $\ell(\cdot, y)$ and $\ell(x, \cdot)$ satisfy the RSI with $\mu = 1/4$.

*Proof.* Proof. We want to show that

$$ (2x + \sin(2x)(3 + a \sin^2(y))/2)(x - x^*) \geq \mu(x - x^*)^2. \quad (11) $$

First note that for any $y$, $\ell(x, y)$ has a unique minimizer $x^* = 0$. Further, for $x = 0$, (11) holds for any $\mu$. Hence, assume $x \neq 0$. Then (11) is satisfied if

$$ 2 + \frac{\sin(2x)}{2x}(3 + a \sin^2(y)) \geq \mu, $$

Using the fact that $\sin(x) \in [-1, 1]$ and that the minimum of $\sin(x)/x$ is attained at $\{-\frac{3}{2}\pi, \frac{3}{2}\pi\}$ (see, e.g., Klén et al. [2010]), $\sin(2x)/(2x) \geq -0.22$. Therefore, since $3 \leq (3 + a \sin^2(y)) \leq 4$ for any $a$ and $y$, we have

$$ \frac{\sin(2x)}{2x}(3 + a \sin^2(y)) \geq -0.88, \quad \text{so} \quad 2 + \frac{\sin(2x)}{2x}(3 + a \sin^2(y)) \geq 1.12, $$

proving the claim.

**Proposition B.2.** The game defined in (1) is not convex and a fortiori not strongly convex.

*Proof.* Proof. The second derivative with respect to $x$ evaluated at $x, y$ is $2 + \cos(2x)(3 + a \sin^2(y))$. which can be negative, e.g. at $x = \pi/2$ and $y = 0$.

**Proposition B.3.** Contraction is satisfied for problem (1) with GDA with $C = (1 - \mu_{\text{RSI}}^2/L^2)^{1/2}$ where $\mu_{\text{RSI}}$ denotes the RSI parameter.

*Proof.* Proof. As for the strongly-convex case:

$$ \|x_2 - x^*\|^2 = \|x - x^* - \eta \nabla_x \ell(x, y)\|^2 = \|x - x^*\|^2 - 2\eta \langle x - x^*, \nabla_x \ell(x, y) \rangle + \eta^2 \|\nabla_x \ell(x, y)\|^2. $$

By RSI

$$ \langle x - x^*, \nabla_x \ell(x, y) \rangle + \langle y - y^*, -\nabla_y \ell(x, y) \rangle \geq \mu_{\text{RSI}} \|x - x^*\|^2 + \mu_{\text{RSI}} \|y - y^*\|^2, $$

thus,

$$ \|x_2 - x^*\|^2 + \|y_2 - y^*\|^2 \leq (1 - 2\eta\mu_{\text{RSI}})(\|x_1 - x^*\|^2 + \|y_1 - y^*\|^2) + \eta^2(\|\nabla_x \ell(x, y)\|^2 + \|\nabla_y \ell(x, y)\|^2). $$

The result follows as before via L-smoothness.
**Lemma C.1.** For any sequence of comparators ($C_t$), the meta-algorithm with adaptive learning rate and losses as defined in (7) and (8) ensures that

$$ \sum_{t=1}^{T} \left( \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 \right) \le 2 \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+ + (\hat{G}_T + \mu \hat{D}_T) \sum_{t=1}^{T} \|Z_t^{(k)} - C_t\|. $$

*Proof.* We start by applying the standard exponential weight analysis, which relates the cumulative mix loss to the minimal cumulative loss among the K base algorithms. (We refer the reader to Lemma 1, paragraph 3 and Lemma 2 in De Rooij et al. [2014] for a proof.) This yields

$$ \sum_{t=1}^{T} m_s \le \min_{i \in [K]} \sum_{t=1}^{T} \ell_{s,i} + \frac{\log K}{\lambda_T}. $$

Therefore, using the fact that the $u \le u_+$ for any $u \in \mathbb{R}$,

$$ \sum_{t=1}^{T} \bar{\ell}_t = \sum_{t=1}^{T} (\bar{\ell}_t - m_t) + \sum_{t=1}^{T} m_t \le \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+ + \min_{i \in [K]} \sum_{t=1}^{T} \ell_{t,i} + \frac{\log K}{\lambda_T}. \quad (13) $$

Next we use that by definition $(\log K)/\lambda_T = \sum_{t=1}^{T-1} (\ell_t - m_t)_+$. Plugging this into (13) and rearranging the terms gives

$$ \sum_{t=1}^{T} \bar{\ell}_t - \min_{i \in [K]} \sum_{t=1}^{T} \ell_{t,i} \le 2 \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+. $$

Replacing $\bar{\ell}_t$ and $\ell_{t,i}$ by their definitions in terms of $F_t$ and $Z_t^{(i)}$, we obtain

$$ \sum_{t=1}^{T} \langle F_t(Z_t), Z_t \rangle - \min_{i \in [K]} \sum_{t=1}^{T} \left( \langle F_t(Z_t), Z_t^{(i)} \rangle + \frac{\mu}{2} \|Z_t^{(i)} - Z_t\|^2 \right) \le 2 \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+. $$

Let us now rearrange the terms and subtract $\sum_{t=1}^{T} \frac{\mu}{2} \|Z_t - C_t\|^2$ on both sides to obtain

$$ \begin{aligned} & \sum_{t=1}^{T} \left( \langle F_t(Z_t), Z_t - C_t \rangle - \frac{\mu}{2} \|Z_t - C_t\|^2 \right) \\ & \le \min_{i \in [K]} \sum_{t=1}^{T} \left( \langle F_t(Z_t), Z_t^{(i)} - C_t \rangle + \frac{\mu}{2} \|Z_t^{(i)} - Z_t\|^2 - \frac{\mu}{2} \|Z_t - C_t\|^2 \right) + 2 \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+ \\ & \le \sum_{t=1}^{T} \left( \langle F_t(Z_t), Z_t^{(k)} - C_t \rangle + \frac{\mu}{2} \|Z_t^{(k)} - Z_t\|^2 - \frac{\mu}{2} \|Z_t - C_t\|^2 \right) + 2 \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+. \end{aligned} $$

Finally, for any $t$, the loss of the k-th base algorithm can be controlled by the distance to the comparators:

$$ \begin{aligned} & \langle F_t(Z_t), Z_t^{(k)} - C_t \rangle + \frac{\mu}{2} \|Z_t^{(k)} - Z_t\|^2 - \frac{\mu}{2} \|Z_t - C_t\|^2 \\ & = \langle F_t(Z_t), Z_t^{(k)} - C_t \rangle + \frac{\mu}{2} \langle Z_t^{(k)} + C_t - 2Z_t, Z_t^{(k)} - C_t \rangle \le (\hat{G}_T + \mu \hat{D}_T) \|Z_t^{(k)} - C_t\|, \end{aligned} $$

where we used the Cauchy-Schwarz inequality to conclude. Summing over $t$ gives the claim. ■

We now proceed to bound the cumulative mixability gap separately. The proof is essentially an application of the mixability of the square loss with respect to the mean, which is equivalent to its exp-concavity. Recall a function $f$ is said to be $\alpha$-exp-concave if $\exp(-\alpha f)$ is convex. For completeness, we provide a self-contained proof of the following result, with no explicit call to mixability.

**Lemma C.2.** Under the assumptions of Lemma C.1,

$$ \sum_{t=1}^{T} (\bar{\ell}_t - m_t)_+ \le \frac{2(\mu \hat{D}_T + \hat{G}_T)^2}{\mu} \log K + \hat{G}_T \hat{D}_T + \frac{\mu}{2} \hat{D}_T^2. $$