Sampling from the Gibbs Distribution in Congestion Games

1. Introduction

Congestion games constitute a rich class of games that have been studied extensively since their introduction by Rosenthal [55]. An (unweighted) congestion game $\mathrm{\Gamma }=(N,E,{({\mathcal{S}}_i)}_{i\in N},{(c_e)}_{e\in E})$ consists of a set of players $N=\{1,\dots ,n\}$ and a set of resources $E=\{1,\dots ,m\}$. Every player i has a strategy set ${\mathcal{S}}_i\subseteq 2^E$, in which each strategy is a subset of resources. Furthermore, every resource $e\in E$ is equipped with a cost function $c_e:{\mathbb{R}}_{\ge 0}\to \mathbb{R}$ that we assume to be nonnegative and nondecreasing. The goal of a player is to choose a strategy that minimizes the player’s total cost $C_i(s)={\sum }_{e\in s_i}{c}_e({\ell }_e(s))$, where ${\ell }_e(s)$ is the number of players using resource e in profile $s\in {\times }_i{\mathcal{S}}_i=\mathcal{S}$. A well-known example is the class of symmetric network congestion games, in which we are given a directed graph $G=(V,E)$ with origin $o\in V$ and destination $d\in V$. The common strategy set of all players is given by the set of all o, d-paths in G.

Rosenthal [55] proves that congestion games are (exact) potential games. He shows that the function $\mathrm{\Phi }:{\times }_i{\mathcal{S}}_i\to \mathbb{R}$ given by

satisfies for every $s\in {\times }_i{\mathcal{S}}_i$ and $s_i^{\prime }\in {\mathcal{S}}_i$ the equality

Here, $C_i(s_i^{\prime },s_{-i})$ is used to denote the cost of player i in the strategy profile in which i chooses $s_i^{\prime }$, and all other players choose their strategy in s. The function $\mathrm{\Phi }$ is often referred to as Rosenthal’s potential. The main implication of (1) is the existence of a pure Nash equilibrium (PNE): a strategy profile in which no player can deviate to another strategy and obtain an improved cost (Rosenthal [55]). This follows directly from the observation that better (or best) response dynamics converge to a PNE in a finite number of steps. Better response dynamics is defined as the procedure by which, in every step, precisely one player deviates to another strategy that yields an improved cost (until a pure Nash equilibrium is reached). For best response dynamics, the deviating player always deviates to a strategy that yields the greatest possible improvement in cost.

In the last two decades, the algorithmic aspects of (pure) Nash equilibria are studied extensively in both general and special classes of congestion games. Two of the most prominent questions concerning pure Nash equilibria are the following.

1. Does (natural) player dynamics, such as better or best response dynamics, converge to a PNE in polynomial time?

2. If not, can one compute a PNE in polynomial time by other means?

Ieong et al. [38] consider the convergence time of better and best response dynamics for singleton congestion games. Here, every strategy of a player consists of a single resource. It is shown that better response dynamics converge in at most $O(n^{2}m)$ steps, and there exist instances in which convergence of better response dynamics might take $\mathrm{\Omega }(n^{\frac{3}{2}}m)$ steps (Ieong et al. [38]). The result of Ieong et al. [38] is generalized by Ackermann et al. [2] to matroid congestion games, in which the strategy set of every player is the collection of bases of a matroid on the ground set E (see Section 2.3 for a formal definition). They show that best response dynamics converge in at most $O(n^{2}mr)$ steps, where r is the maximum rank of the matroids that form the strategy sets of the players. They also show that matroids are, in a sense, the maximal structure that allow for polynomial time convergence of best response dynamics. Furthermore, Ackermann et al. [2] show best response dynamics are not guaranteed to converge in symmetric network congestion games, in which the common strategy set of every player is the collection of paths in a given directed graph. Best response dynamics are also studied in more general congestion games, such as with player-specific cost functions (Ackermann and Röglin [1]) or with weighted players (Even-Dar et al. [24], Fanelli and Moscardelli [26]). Approximate versions of response dynamics are also studied; see, for example, Chien and Sinclair [19] and Skopalik and Vöcking [59].

The complexity of computing a pure Nash equilibrium is also well-understood. Fabrikant et al. [25] show that computing a pure Nash equilibrium is complete for the complexity class polynomial local search. Hardness of computing a PNE is also shown for various special cases of congestion games; see, for example, Ackermann et al. [2] and Del Pia et al. [21], as well as for the computation of approximate equilibria (Skopalik and Vöcking [59]). On the positive side, Fabrikant et al. [25] show that a PNE can be computed efficiently in symmetric network congestion games despite the fact that best response dynamics are not guaranteed to converge in these games (Ackermann et al. [2]). Del Pia et al. [21] generalize this result to totally unimodular congestion games; see also Kleer and Schäfer [43] for a more general framework.

In general, player dynamics roughly comes in two flavors: one deviates to another strategy profile according to either a deterministic rule or a probabilistic one. A well-known example of the latter case is noisy (randomized) best response dynamics, which has received a lot of attention in practice; see, for example, Camerer [18] and Mäs and Nax [47], but seems hard to analyze from a theoretical perspective. Here, instead of making a deviation to another strategy according to a deterministic rule, a player chooses a strategy from the player’s set according to a probability distribution that usually puts relatively more weight on strategies that result in a lower cost. Randomized dynamics can be studied from two perspectives: as either a randomized alternative for deterministic dynamics converging to a pure Nash equilibrium or as a dynamical system on its own.

One well-known example of player dynamics that can be studied as a dynamical system and which is the topic of this paper is logit dynamics. It has received a lot of attention in various communities, such as evolutionary game theory (e.g., Sandholm [56]) and experimental economics (e.g., Camerer [18]). The procedure was introduced by Blume [14] as a form of randomized game dynamics in which players update their strategy according to a logit update rule (McFadden [48]). The logit dynamics for congestion games can be formulated as follows. For a given strategy profile s and fixed rationality level (or inverse temperature in the physics literature) parameter $T\ge 0$, first choose a player $i\in N$ uniformly at random and then have player i choose a strategy $s_i^{\prime }\in {\mathcal{S}}_i$ with probability

Note that the denominator in (2) is a normalizing constant (often called the partition function).

Remark 1.

Equivalently, one can replace the $\mathrm{\Phi }(s_i^{\prime },s_{-i})$ by $C_i(s_i^{\prime },s_{-i})$ and, similarly, in the normalizing constant. We use $\mathrm{\Phi }$ as this is more convenient for our purposes. The equivalence follows from (1).

The rationality level $T\ge 0$ is used to model the amount of noise players believe there to be in the system (Auletta et al. [12]). When $T\to \infty$, players effectively only assign positive probability to best responses, whereas when $T\to 0$, the distribution in (2) approaches the uniform distribution over ${\mathcal{S}}_i$. Also note that the dynamics indeed puts relatively more probability mass on strategies that give a greater improvement in cost.

Logit dynamics gives rise to an ergodic, time-reversible Markov chain on the set $\mathcal{S}$ of all strategy profiles that has the Gibbs distribution π given by

for $s\in \mathcal{S}$ as its unique stationary distribution. This simply means that, if one runs the logit dynamics for a sufficiently long time, the distribution over $\mathcal{S}$ converges to the Gibbs distribution. The time it takes to converge to the stationary distribution is called the mixing time of the Markov chain.

Auletta et al. [12] interpret the Gibbs distribution as a dynamic equilibrium concept, which they dubbed the logit equilibrium; see also Ferraioli [28]. This concept is well-defined for general finite games; see, for example, Auletta et al. [12]. The goal of this work is to study algorithmic aspects of the logit equilibrium/Gibbs distribution in congestion games. Auletta et al. [12] give a tight bound for the mixing time of the logit dynamics in general potential games, which can be exponential in the quantity $T{\mathrm{\Phi }}_{\max }$ (as in Example 1). They also study graphical coordination games, for which they identify various types of behavior of the logit dynamics; we refer to Auletta et al. [12] for more details. Because of the general slow mixing, Auletta et al. [10] also study the concept of metastability, which means a Markov chain stays close to some distribution on a timescale shorter than the mixing time. More work on the logit dynamics in potential games by these authors can also be found in Auletta et al. [9, 11]. We also refer to the survey article of Ferraioli [28] and references therein. There are also various results addressing the inefficiency of “long-term” equilibria in the context of logit dynamics; see, for example, the works of Asadpour and Saberi [8], Mamageishvili and Penna [45], and Penna [54].

A special case of potential games for which logit-like dynamics is studied extensively is the Glauber dynamics for the Ising model, which, in game-theoretical terms, can be seen as logit dynamics for max-2-cut games that are described in Section 1.1. Whether the logit dynamics are rapidly mixing in this case depends on the parameter $T\ge 0$ and graph topology G; see, for example, the work of Levin et al. [44] and references therein. Jerrum and Sinclair [39] show that, nevertheless, there exists a polynomial time algorithm to sample from the Gibbs distribution for any rationality level $T\ge 0$ and any graph topology. This work is partially motivated by studying a similar question for special cases of congestion games as explained later.

As said before, in general, the logit dynamics converge slowly to the Gibbs distribution (Auletta et al. [12]); in particular, the number of steps needed might be $\mathrm{\Omega }(e^{T{\mathrm{\Phi }}_{\max }})$, where ${\mathrm{\Phi }}_{\max }$ is the maximum value attained by Rosenthal’s potential. We next give a simple example illustrating this fact.

Example 1.

Consider the congestion game with players $N=\{1,2\}$ and two resources $E=\{a,b\}$. Both resources $e\in \{a,b\}$ can be used by both players, that is, ${\mathcal{S}}_1={\mathcal{S}}_2=\{\{a\},\{b\}\}$, and have a cost function satisfying $c_e(1)=0$ and $c_e(2)=\varphi$ for some $\varphi \ge 0$. Note that ${\mathrm{\Phi }}_{\max }=\varphi$.

We assume that T is fixed in this example, in particular, independent of $\varphi$. If $\varphi$ is large, the Gibbs distribution assigns weight close to 1/2 to the strategy profiles $(s_1,s_2)\in \{(a,b),(b,a)\}$ and weight close to (1/2)$e^{-T\varphi }\approx 0$ to $(s_1,s_2)\in \{(a,a),(b,b)\}$. We summarize these facts in Figure 1. Now, informally speaking, if we consider the logit dynamics with starting profile (a, b), then in order to reach the profile (b, a), we have to go through either the profile (a, a) or (b, b). The probability of transitioning to one of those profiles in one step of the logit dynamics is $O(e^{-\varphi T})$. As both profiles (a, b) and (b, a) appear with probability close to 1/2 in the Gibbs distribution, the probability of being in either of these profiles after a (large) number of steps of the logit dynamics should be approximately equal. However, escaping from the profile (a, b) happens with probability $O(e^{-\varphi T})$. The inverse $\mathrm{\Omega }(e^{T\varphi })$ of this probability then forms a lower bound on the number of steps of the logit dynamics needed to get close to the stationary Gibbs distribution.

What is causing the slow convergence in Example 1? The problem is that we need to use either the profile (a, a) or (b, b), both having very small probability in the Gibbs distribution, to move from (a, b) to (b, a). However, as it turns out, if one in addition with some probability is allowed to interchange the strategies of players 1 and 2, then the resulting dynamics converges quickly to the Gibbs distribution (being a special case of Theorem 1). Note that this enables the possibility to directly transition between (a, b) and (b, a). This motivates the following question: Is there a (simple) Markov chain on $\mathcal{S}$ that converges rapidly to the Gibbs distribution over $\mathcal{S}$ at any rationality level T?

This question may be interpreted as the natural analogue of looking for other local search procedures converging quickly to a PNE when best/better response dynamics does not have this property.

When the answer to the preceding question is not directly obvious, one can take another step back and first ask whether it is at all possible to efficiently sample from the Gibbs distribution. Informally speaking, can we take “snapshots” (according to the Gibbs distribution) from the system in equilibrium in polynomial time? More formally speaking, does there exist an efficient algorithm to sample (approximately) a strategy profile $s\in \mathcal{S}$ according to the Gibbs distribution over $\mathcal{S}$ at any rationality level T? This question can be interpreted as a dynamic analogue of the second question posed earlier for the computation of pure Nash equilibria. That is, although (deterministic) better/best response dynamics might take a long time to converge to a pure Nash equilibrium, one still wants to know whether a PNE can be computed efficiently by other means. Similarly, although logit (or other natural) dynamics can take a long time to converge to the Gibbs distribution, we may still ask whether, by means of sampling, we can get an impression of what the Gibbs distribution over $\mathcal{S}$ looks like. Note that a positive answer to the first question posed earlier gives a positive answer to the second question as one can simply run the Markov chain sufficiently long in order to generate a sample from the Gibbs distribution.

These questions are made precise in Section 2. We remark that one of the aspects that makes them nontrivial is the fact that we require the questions to hold for any rationality level T, that is, T is considered part of the input. For example, in Example 1, we could set $T=\mathrm{\Theta }(1/\varphi )$ to circumvent the problem arising there.

1.1. Our Contributions

In this work, we address the questions posed for various special cases of congestion games.

1.1.1. Extension Parallel Congestion Games.

The class of extension parallel (EP) congestion games is a well-studied special case of symmetric network congestion games; see, for example, Fotakis [29], Fujishige et al. [30], and Holzman and Law-Yone [37]. Here, the common strategy set of all players is given by the set of o, d-paths $\mathcal{P}$ of an extension parallel graph; see Section 4 for a definition and example. Our main result is that there is a simple Markov chain converging quickly to the Gibbs distribution over $\mathcal{S}$, also implying that we can sample approximately from the Gibbs distribution. We show that, if one, in addition to the logit dynamics transitions, is allowed to randomly interchange the strategies of two players (akin to the explanation given after Example 1), the resulting Markov chain converges quickly to the Gibbs distribution. We call this the relaxed logit dynamics; see Section 4.2 for a formal definition. Note that Theorem 1 gives a doubly exponential improvement with respect to the dependence on $T{\mathrm{\Phi }}_{\max }$ compared with the lower bound as given in Example 1.

Theorem 1

(Informal). The relaxed logit dynamics for EP congestion games, at rationality level T, converges to a distribution ϵ-close to the Gibbs distribution in at most

$$n^3\left(\log n+\log \log |\mathcal{P}|+\log \left(\frac{2T{\mathrm{\Phi }}_{\max }}{{\mathrm{ϵ}}^2}\right)\right)$$

steps, where n is the number of players, $\mathcal{P}$ the number of paths in the EP graph, and ${\mathrm{\Phi }}_{\max }$ the maximum value attained by Rosenthal’s potential.¹

The notion of “ϵ-close” refers to the fact that the distribution seen after the indicated number of steps differs from the Gibbs distribution at most ϵ in total variation distance (see Section 2.5), a well-known distance measure for comparing probability distributions in Markov chain theory.

In a nutshell, Theorem 1 follows from the fact that, in EP congestion games, Rosenthal’s potential is M-convex as is shown by Fujishige et al. [30]. M-convexity is a property defined in the area of discrete convex analysis (Murota [51]) (see Section 2.3).² The link between M-convexity and sampling is, roughly speaking, established in a series of (breakthrough) papers by Anari et al. [3, 5, 6] through the theory of strongly log-concave (SLC) polynomials. The theory of strongly log-concave (or Lorentzian) polynomials dates back to the work of Gurvits [34] and is further developed by Anariet al. [3] and Brändén and Huh [17]. In particular, Anari et al. [5] give the first polynomial-time algorithm for approximately sampling and counting the number of bases of a given matroid, resolving also an old conjecture by Mihail and Vazirani [50]. In this work, we rely on the sampling result from Anari et al. [5], albeit for relatively simple matroid structures.

Before proving the preceding result, we also give another way of sampling from the Gibbs distribution in Section 4.1 that essentially is a more direct approach than the sampler induced by the Markov chain result given (Theorem 4). The high-level approach used for this more direct sampler is given in Section 3 (and also used in Section 6). Finally, we give an application of our results to the problem of (approximate) uniformly sampling pure Nash equilibria in EP congestion games in Section 4.3, which is, to the best of our knowledge, the first of its kind.

1.1.2. Max-k-Cut Game on a (Unweighted) Complete Graph.

A max-k-cut game (see, e.g., Gourvès and Monnot [33]) is given by an undirected graph $G=(V,E)$ whose nodes V are the players. Every $v\in V$ has strategy set ${\mathcal{S}}_v=\{1,\dots ,k\}$ whose elements are referred to as colors. The utility of a player is the number of neighbors that choose a different color. They can be modeled as a special case of extension parallel congestion games when G is the complete graph. This is explained in Section 5.

As an application of Theorem 1, we show that the relaxed logit dynamics is rapidly mixing at any rationality level T for max-k-cut games when G is the complete graph (Corollary 2). This stands in stark contrast to known results regarding the logit dynamics on complete graphs for the case k = 2, whose mixing time depends heavily on the value of the rationality level T (Levin et al. [44]). We elaborate on this in Section 5.

1.1.3. Capacitated Uniform Congestion Games.

Finally, we study the class of so-called u-capacitated k-uniform congestion games for given $k=(k_1,\dots ,k_n)$ and $u=(u_1,\dots ,u_m)$. In such a game, the strategy set of player $i\in N$ is given by all subsets of E of size k_i. Furthermore, for every $e\in E$, we are given a capacity u_e so that $c_e(x)=\infty$ whenever $x>u_e$.

The motivation for studying these games comes from the class of base-matroid congestion games, in which the strategy set of every player is the set of bases ${\mathcal{B}}_i$ of a given matroid ${\mathcal{M}}_i$ over the ground set of resources E. It is well-known that best response dynamics converge to a PNE in a polynomial number of steps in this class of games (Anari et al. [2]), and so, in particular, a PNE can be computed in polynomial time. Given the base-matroid sampling result of Anari et al. [5], a natural question that comes to mind is if a similar result exists for sampling from the Gibbs distribution in base-matroid congestion games. Here, we give a first result addressing this question. (Note that, in our setting, the strategies of player $i\in N$ are the bases of the k_i-uniform matroid.)

We next explain why there is a need for capacity constraints in our results. A strategy profile s can be seen as a bipartite graph in which the nodes on one side correspond to the players, having degrees k_i, and the nodes on the other side correspond to the resources, having degrees ${\ell }_e(s)$ (the resource load on e in profile s). For a given profile s with resource load profile $\ell (s)$, the bipartite graph is obtained in the natural way: there is an edge between player i and resource e if player i uses resource e in profile s. Very roughly speaking, in order to apply the sampling result in Anari et al. [5], we establish strong log-concavity of a certain polynomial associated to the vectors k and u. For this, we rely on an asymptotic enumeration formula for the number of bipartite graphs with a given degree sequence, which, in terms of congestion games, gives the number of strategy profiles that have a given resource load profile $\ell$. (Asymptotic enumeration of graphs with given degrees is studied extensively in the area of combinatorics.) The formula that we use is only valid for the range of $k=(k_1,\dots ,k_n)$ and resource load profiles $\ell (s)=({\ell }_e(s))$ satisfying the imposed capacity constraints as given in Theorem 2. Our main result is as follows.

Theorem 2

(Informal). There is an (almost) polynomial time algorithm for approximately sampling from the Gibbs distribution in u-capacitated k-uniform congestion games assuming that $1\le k_{\max }{u}_{\max }=o\left(U^{1/4}\right)$ when $n\to \infty$, where $k_{\max }={\max }_i{k}_i,\hspace{0.17em}{u}_{\max }={\max }_j{u}_j$ and $U={\sum }_j{u}_j$.

See Remark 6 for an explanation of what we mean by “almost” polynomial time. The proof of Theorem 2 reveals an interesting connection between M-convexity and asymptotic enumeration formulas that might be of independent interest.

1.2. Technical Approach

In this section, we give an intuitive outline of the approach taken in order to prove the results as summarized in Section 1.1. For the sake of simplicity, we work with symmetric singleton congestion games in which the common strategy set contains every individual resource as a strategy, that is, ${\mathcal{S}}_i=\{\{e_1\},\{e_2\},\dots ,\{e_m\}\}$ for all $i\in N$. These games form a special case of EP and capacitated uniform congestion games. All technical notions not formally defined here can be found in Section 2.

Let us recall Example 1, which is a symmetric singleton congestion game with two resources. The cause of the long convergence time in this example is the fact that transitioning between strategy profiles in which both resources have the same load profile, that is, the number of players using the resource, is difficult. To be precise, in Example 1, both profiles (a, b) and (b, a) have the same load profile, but they correspond to different strategy profiles.

In order to still be able to sample from the Gibbs distribution, we take two approaches:

1. We first sample a load profile according to the correct probability and then sample a strategy profile corresponding to that load profile (Algorithm 1 in Section 3).

2. We augment the logit dynamics with additional transition possibilities, making it easier to transition between strategy profiles with the same load profile.

In the first approach, the “correct probability” with which a given load profile should be sampled is the sum of all the probabilities assigned to strategy profiles with the given profile. The resulting distribution is called the Gibbs distribution induced on the load profiles. For the applications in Sections 4–6, sampling a strategy profile with a given load profile can be done by either generating a random permutation of the players or relying on known sampling results in the literature. In particular, in the case of capacitated uniform congestion games, sampling a strategy profile for a given load profile corresponds to sampling a bipartite graph with a given degree sequence.

The main technical challenge lies in sampling a load profile with the correct probability. In the case of symmetric singleton congestion games, the collection of all load profiles can be modeled by the set $R=\{x\in {\mathbb{Z}}_{\ge 0}^m:x_1+x_2+\cdots +x_m=n\}$, where x_e denotes the load on resource $e\in E$. The set R is a special case of a discrete polymatroid. In order to sample a (resource) load vector x from R with the correct probability, we consider a Markov chain that can transition between load vectors x and y for which ${\sum }_i|x_i-y_i|=2$; that is, x and y differ in one unit of load on precisely two resources (but have the same total load). The transition probabilities are chosen in such a way that the stationary distribution of the Markov chain is precisely the Gibbs distribution induced on the load profiles.

The Markov chain that we consider over the set R is a “polymatroid version” of the base-exchange Markov chain studied by Anari et al. [5] for sampling a base of a given matroid. The state space of this base-exchange Markov chain is the set of all bases of a given matroid. It randomly exchanges some element of the current base with another element not in the base, resulting in another base of the matroid. This idea relies on the well-known base-exchange property of matroids (Schrijver [57]). In their breakthrough work Anari et al. [5] show this Markov chain to be rapidly mixing when the stationary distribution is strongly log-concave. Strong log-concavity is a property of a polynomial associated with the stationary distribution of a Markov chain (defined for both discrete polymatroids and matroids). We apply the result from Anari et al. [5] by giving a reduction from sampling a base of a polymatroid, that is, a load vector from R, to the problem of sampling a base of a matroid using a well-known construction of Helgason [35]. For this reduction, we rely on some properties that preserve strong log-concavity (Brändén and Huh [17]). Whereas this reduction works for arbitrary discrete polymatroids, in our applications, we use relatively simple polymatroid structures. What is left then is to show that the Gibbs distribution induced on the load profiles indeed satisfies the notion of strong log-concavity. Based on results in Brändén and Huh [17], it turns out that, for symmetric singleton congestion games, this is the case because of the fact that Rosenthal’s potential is a separable convex function.

In order to obtain our sampler for EP congestion games (Theorem 4), we use a similar approach based on load profiles, but instead, we rely on the M-convexity of Rosenthal’s potential for such games (Fujishige et al. [30]). For Theorem 2, we also use a similar approach.

In order to show that the logit dynamics converge quickly to the Gibbs distribution when we additionally allow the interchanging of the strategy of two players, we use a Markov chain decomposition argument. This approach results in the proofs of Theorem 5 and Corollary 2. The idea of Markov chain decomposition is to divide the state space into smaller sets. If one can prove that the Markov chain converges quickly to its (induced) stationary distribution on the smaller sets and that it is easy to move between these smaller sets, then the original chain converges quickly to its stationary distribution as well. We elaborate more on this approach before giving the proof of Theorem 5.

1.3. Discussion and Further Related Work

To the best of our knowledge, ours are the first (polynomial-time) sampling results for the Gibbs distribution in congestion games beyond the well-studied case of max-cut games on general graphs, also known as the Glauber dynamics for the Ising model. Given the extensive attention that logit dynamics has received in various communities as well as the topic of approximate sampling in theoretical computer science, we believe this to be an interesting line of work at the intersection of algorithmic game theory, combinatorics, and approximate sampling to pursue further.

In particular, for special cases of congestion games with a positive answer to questions 1 and 2 posed in the introduction, do there also exist positive answers for their dynamic analogues? As a concrete open question, we ask whether it is always possible to efficiently sample from the Gibbs distribution in general base matroid congestion games (Ackermann et al. [2]). If true, this provides an interesting (qualitative) game-theoretical generalization of the sampling result of Anari et al. [5].

We conclude this section with some more background and related work concerning the result of Anari et al. [5]. A long-standing open question (Mihail and Vazirani [50]) in the area of approximate sampling and counting asks if there exists a rapidly mixing Markov chain to sample approximate uniformly at random a base of a given matroid (which, in turn, also yields a polynomial-time approximation scheme for counting the number of bases of a given matroid).³ Anari et al. [5] show that the base-exchange Markov chain is indeed rapidly mixing. In fact, they prove the more general result that this chain is rapidly mixing, using appropriate transition probabilities, for any strongly log-concave stationary distribution; see Section 2.4. An (up to constant factors) tight mixing time bound for the base-exchange Markov chain is later given by Anari et al. [6]. Before Anari et al. [5], rapid mixing of the base-exchange Markov chain (with uniform stationary distribution) was only known for special cases. In particular, Feder and Mihail [27] prove rapid mixing for balanced matroids, using Sinclair’s [58] multicommodity flow method. Examples of balanced matroids are the graphic and regular matroids.

2. Preliminaries

In this section, we give all the necessary preliminaries regarding congestion games, strongly log-concave polynomials, and the relevant Markov chain notions and results. We start with some general notation.

All logarithms in this work have Euler’s number e as their base unless specified otherwise. For $k\in {\mathbb{Z}}_{>0}$, we write $[k]=\{1,\dots ,k\}$. For two vectors $x,y\in {\mathbb{Z}}^n$, we write $x\le y$ if $x_i\le y_i$ for $i=1,\dots ,n$, and x < y if strict inequality holds for at least one i. Furthermore, with $|x|={\sum }_{i=1}^n|x_i|$, we denote the modulus of x. We use ${(e_i)}_{i=1,\dots ,n}$ to denote the standard basis of ${\mathbb{R}}^n$, that is, $e_i(j)=1$ if j = I and $e_i(j)=0$ otherwise. For sets A, B, and C, we use the notation $A-B+C$ to denote the set $(A\backslash B)\cup C$.

2.1. Congestion Games

A capacitated congestion game Γ is given by a tuple $(N,E,{({\mathcal{S}}_i)}_{i\in N},{(c_e)}_{e\in E},{(u_e)}_{e\in E})$, where $N=[n]$ is a finite set of players, $E=[m]$ a finite set of resources (or facilities), ${\mathcal{S}}_i\subseteq 2^E$ the set of strategies of player $i\in N$, and $c_e:{\mathbb{Z}}_{\ge 0}\to \mathbb{Q}$ the cost function of resource $e\in E$ that satisfies $c_e(x)=W$ whenever $x>u_e$ for $e\in E$ with W a sufficiently large number. Unless stated otherwise, the cost functions are assumed to be nonnegative and nondecreasing. Finally, u_e is a nonnegative integer modeling the capacity on resource $e\in E$. If u_e = n for every resource $e\in E$, we simply call Γ a congestion game. (Extension parallel and capacitated uniform congestion games are defined in Sections 4 and 6, respectively.) For a strategy profile $s=(s_1,\dots ,s_n)\in {\times }_i{\mathcal{S}}_i=\mathcal{S}$, we define ${\ell }_e(s)$ to be the number of players using resource e, that is, ${\ell }_e(s)=|\{i\in N:e\in s_i\}|$. A game is called symmetric if ${\mathcal{S}}_i={\mathcal{S}}_j$ for all $i,j\in [n]$. We then write $\mathcal{P}$ to denote the common strategy set of all players.

We call $\ell (s)={({\ell }_e(s))}_{e\in E}$ the resource load profile corresponding to strategy profile s. We say that a strategy $s\in \mathcal{S}$ is feasible if ${\ell }_e(s)\le u_e$ for every $e\in E$ and write ${\mathcal{S}}_f$ to denote the set of all feasible strategy profiles. We only consider games in which the set of feasible strategy profiles is nonempty. More generally, we say that $y\in {\mathbb{N}}^m$ is a (feasible) resource load profile for $(N,E,{({\mathcal{S}}_i)}_{i\in N},{(u_e)}_{e\in E})$ if there is some (feasible) strategy profile s such that $y=\ell (s)$. We write $\mathcal{S}(y)$ for the set of strategy profiles $s\in {\times }_i{\mathcal{S}}_i$ whose resource load profile is y.

Similarly, for symmetric congestion games, we define the notion of a strategy load profile that models how many players are using a strategy $p\in \mathcal{P}$ in a given strategy profile $s\in {\mathcal{P}}^n$. More precisely, given a strategy profile $s\in {\mathcal{P}}^n$, we define $z_p(s)=|\{i\in N:s_i=p\}|$ as the number of players choosing strategy $p\in \mathcal{P}$ in strategy profile s. The vector $z(s)={(z_p(s))}_{p\in \mathcal{P}}$ is called the strategy load profile of s. Similarly as for resource load profiles, we define (with a slight abuse of notation) $\mathcal{S}(x)$ to be the set of strategy profiles s for which $x=z(s)$.

The cost of player $i\in N$ under a strategy profile $s=(s_1,\dots ,s_n)\in {\times }_i{\mathcal{S}}_i$ is given by $C_i(s)={\sum }_{e\in s_i}{c}_e({\ell }_e(s)).$ A strategy profile $s\in \mathcal{S}$ is called a (pure) Nash equilibrium if, for every $i\in N$ and every $s_i^{\prime }\in {\mathcal{S}}_i$, it holds that $C_i(s)\le C_i(s_i^{\prime },s_{-i})$, where $(s_i^{\prime },s_{-i})$ denotes the strategy profile in which player i plays $s_i^{\prime }$ and every other player $j\ne i$ plays s_j. We write $\text{NE}(\mathrm{\Gamma })$ to denote the set of all pure Nash equilibria of Γ.

We say that $\mathrm{\Phi }:{\times }_i{\mathcal{S}}_i\to \mathbb{R}$ is an exact potential function for a congestion game Γ if, for every strategy profile $s\in {\times }_i{\mathcal{S}}_i$, for every player $i\in N$, and every unilateral deviation $s_i^{\prime }\in {\mathcal{S}}_i$ of i it holds that $\mathrm{\Phi }(s)-\mathrm{\Phi }(s_{-i},s_i^{\prime })=C_i(s)-C_i(s_{-i},s_i^{\prime }).$ Rosenthal [55] shows that

is an exact potential function for any congestion game. Subsequently, we refer to this potential function as Rosenthal’s potential.

A function $\varphi :\{0,\dots ,n\}\to \mathbb{R}$ is called convex if $\varphi (i)-\varphi (i-1)\le \varphi (i+1)-\varphi (i)$ for all $i=1,\dots ,n-1$. A function $\psi :{\{0,\dots ,n\}}^m\to \mathbb{R}$ is called separable convex if it is of the form $(x_1,\dots ,x_m)\mapsto {\sum }_{j=1}^m{\psi }_j(x_j)$, where the ψ_j, given by $x_j\mapsto {\psi }_j(x_j)$, are convex. We say that $\varphi$ is concave if $-\varphi$ is convex, and similarly, ψ is separable concave if $-\psi$ is separable convex. A simple but important observation that we use in this work is the fact that Rosenthal’s potential is a separable convex function when seen as a function from resource load profiles to the reals. That is, the function $\tilde{\mathrm{\Phi }}:{\mathbb{Z}}_{\ge 0}^m\to \mathbb{R}$ given by

for $\beta \in {\mathbb{Z}}_{\ge 0}^m$ is separable convex.

Proposition 1.

If the cost functions ${(c_e)}_{e\in E}$ are nondecreasing, then Rosenthal’s potential $\tilde{\mathrm{\Phi }}$ is a separable convex function.

2.2. Gibbs Distribution and Logit Dynamics

The Gibbs distribution $\pi :\mathcal{S}\to {\mathbb{R}}_{\ge 0}$ over the strategy profiles of a congestion game Γ is given by

where $\mathrm{\Phi }$ is Rosenthal’s potential, $T\ge 0$ rationality level, and Z is the normalizing constant (or partition function)

The logit dynamics Markov chain (formal Markov chain definitions are given in Section 2.5) with current state $s\in \mathcal{S}$ proceeds as

• Select a player $i\in N$ uniformly at random.

• For $s_i^{\prime }\in {\mathcal{S}}_i$, transition to $(s_i^{\prime },s_{-i})$ with probability $e^{-T\mathrm{\Phi }(s_i^{\prime },s_{-i})}/Z_i^{\prime }$ with normalizing constant

  $$Z_i^{\prime }={\displaystyle \sum _{r\in {\mathcal{S}}_i}{e}^{-T\mathrm{\Phi }(r,s_{-i})}}.$$

It is a standard fact, which can be shown by using (1), that this Markov chain is reversible with respect to the Gibbs distribution.

2.3. Matroids and M-Concavity

Let $E=[n]$ be a finite set called the ground set and $\mathcal{I}\subseteq 2^E=\{X:X\subseteq E\}$ a collection of subsets of E (called independent sets). The pair $\mathcal{N}=(E,\mathcal{I})$ is a matroid if (i) $\varnothing \in \mathcal{I}$; (ii) $A\in \mathcal{I}$ and $B\subseteq A$, then $B\in \mathcal{I}$; (iii) $A,B\in \mathcal{I}$ and $|A|>|B|$, then there exists an $a\in A\backslash B$ such that $B+a\in \mathcal{I}$. An independent set $B\in \mathcal{I}$ of maximum size is called a basis. We use $\mathcal{B}$ to denote the set of all bases of $\mathcal{N}$. The set of bases $\mathcal{B}$ satisfies the so-called base-exchange property: if $B,B^{\prime }\in \mathcal{B}$ and $e\in B\backslash B^{\prime }$, then there exists an $e^{\prime }\in B^{\prime }\backslash B$ such that $B+e^{\prime }-e\in B$. It satisfies the strong base-exchange property if both $B+e^{\prime }-e,B^{\prime }-e^{\prime }+e\in B$. The rank of a matroid is the common cardinality r of all bases in $\mathcal{B}$. The $\ell$-truncation ${\mathcal{M}}_{\ell }=(E,{\mathcal{I}}_{\ell })$ of a matroid $\mathcal{M}$ is the matroid with $A\in {\mathcal{I}}_{\ell }$ if and only if $A\in \mathcal{I}$ and $|A|\le \ell$. The partition matroid is given by a disjoint partition $E=E_1\cup \dots \cup E_q$ of the ground set E and upper bounds u_i for $i=1,\dots ,q$. A set $A\subseteq E$ is independent if and only if $|A\cap E_i|\le u_i$ for all $i=1,\dots ,q$. The k-uniform matroid is the matroid in which $A\subseteq E$ is independent if and only if $|A|\le k$.

A discrete polymatroid is a finite set of vectors $R\subset {\mathbb{Z}}_{\ge 0}^n$ with the properties that (i) $\mathbf{0}\in R$; (ii) if $y\in R$ and $x\le y$, then $x\in R$; and (iii) if $x,y\in R$ with $|y|>|x|$, then there is a vector $w\in R$ such that $x<w<\max \{x,y\}$ (in which the maximum is taken coordinate-wise). The set of bases ${\mathcal{B}}_R$ is given by all maximal vectors in R that have a common modulus r. A polymatroid satisfies the base-exchange property; if $x,y\in {\mathcal{B}}_R$ and x_i > y_i, then there exists an index j with y_j > x_j and $x-e_i+e_j\in {\mathcal{B}}_R$.

Finally, as a generalization of discrete polymatroids, we describe the notion of M-convexity for functions Murota [51, 53]. As we mostly work with its negated counterpart of M-concavity, we describe this first. Let $\nu :{\mathbb{Z}}_{\ge 0}^n\to \mathbb{R}\cup \{-\infty \}$ be a function. The effective domain of ν is given by

The function ν is called $M^{\mathrm{♯}}$-concave if it satisfies the (symmetric) exchange property: for any $\alpha ,\beta \in \text{dom}(\nu )$ and any $i\in [n]$ satisfying ${\alpha }_i>{\beta }_i$, there exists a $j\in [n]$ such that ${\alpha }_j<{\beta }_j$ and

It is well-known that a separable concave function is $M^{\mathrm{♯}}$-concave (Murota [52]). The function ν is called M-concave if it is $M^{\mathrm{♯}}$-concave and, in addition, there is an $r\in {\mathbb{Z}}_{\ge 0}$ such that $\text{dom}(\nu )\subseteq \{\alpha :{\sum }_i{\alpha }_i=r\}$. A function $\nu :{\mathbb{Z}}_{\ge 0}^n\to \mathbb{R}\cup \{\infty \}$ is called M-convex if $-\nu$ is M-concave.

2.4. Strongly Log-Concave Polynomials

We consider polynomials $p\in \mathbb{R}[x_1,\dots ,x_n]$ with nonnegative real coefficients. For a vector $\beta =({\beta }_1,\dots ,{\beta }_n)\in {\mathbb{Z}}_{\ge 0}^n$, we write

to denote the partial differential operator that differentiates a function β_i times with respect to x_i for $i=1,\dots ,n$. For $\alpha \in {\mathbb{Z}}_{\ge 0}^n$, we write $x^{\alpha }$ to denote ${\prod }_{i=1}^n{x}_i^{{\alpha }_i}$. Furthermore, we write $\alpha !={\prod }_i{\alpha }_i!$, and for $\alpha ,\kappa \in {\mathbb{Z}}_{\ge 0}^n$ with ${\alpha }_i\le {\kappa }_i$ for all i, we write

For a constant $t\in {\mathbb{Z}}_{\ge 0}$ with $t\ge {\max }_i{\alpha }_i$, we write $\left(\begin{array}{c}t\\ \alpha \end{array}\right)={\prod }_{i=1}^n\left(\begin{array}{c}t\\ {\alpha }_i\end{array}\right)$. Let $\kappa \in {\mathbb{Z}}_{\ge 0}^n$ and $K={\times }_i\{0,\dots ,{\kappa }_i\}$. Let $w:K\to {\mathbb{R}}_{\ge 0}$ be a weight function. The generating polynomial of w is given by

The support of $g_{\kappa }$ is the set $\text{supp}(g_{\kappa })=\{\alpha \in K:w(\alpha )>0\}$. The generating polynomial g is called d-homogeneous if $|\alpha |={\sum }_i{\alpha }_i=d$ for all $\alpha \in {\{0,\dots ,k\}}^m$ with $w(\alpha )>0$. It is called multiaffine if every variable x_i appears with at most multiplicity one in every monomial of p. For example, $q(x_1,x_2)=x_1{x}_2$ is multiaffine, but $r(x_1,x_2)=x_1^2+x_1{x}_2$ is not as the multiplicity of x₁ in the first monomial is two. Finally, the elementary symmetric polynomial of degree d, for $\kappa =(1,1,\dots ,1)$, is given by

Definition 1

(Strong Log-Concavity; Gurvits [34]). A polynomial $p\in \mathbb{R}[x_1,\dots ,x_n]$ with nonnegative coefficients is called log-concave on a subset $S\subseteq {\mathbb{R}}_{\ge 0}^n$ if its Hessian ${\mathrm{▽}}^2\log (p)$ is negative semidefinite on S. A polynomial p is called SLC on S if, for any $\beta \in {\mathbb{Z}}_{\ge 0}^n$, we have that ${\partial }^{\beta }p$ is log-concave.

For convenience, the zero polynomial is defined to be strongly log-concave always. It is interesting to note that if a d-homogeneous multiaffine polynomial p is SLC, then the support of p must form the collection of bases of a matroid, and more generally, if a (not multiaffine) homogeneous polynomial is SLC, its support forms an M-convex set (Brändén and Huh [17]). Finally, if the generating polynomial $g_{\kappa }$ is strongly log-concave, then the probability distribution $\pi (\alpha )\propto w(\alpha )$ is called strongly log-concave.

Remark 2.

The definition of strong log-concavity is in this work is not really needed but included for completeness. In our proofs, we essentially only rely on properties of SLC polynomials from the literature that are reviewed as follows. For homogeneous generating polynomials the notion of strong log-concavity is equivalent to that of a polynomial being Lorentzian (Brändén and Huh [17]) or completely log-concave (Anari et al. [4]). These equivalences are shown in Brändén and Huh [16].

We next state all properties of SLC polynomials that are used in this work. First, it is easy to check that the SLC property is preserved under multiplication with a nonnegative scalar, which we state here for sake of reference.

Proposition 2

(Brändén and Huh [17]). If $p\in \mathbb{R}[x_1,\dots ,x_n]$ is SLC and $\gamma \in {\mathbb{R}}_{\ge 0}$, then $\gamma ·p$ is SLC.

We continue with the polarization operator defined by Brändén and Huh [17]. The polarization operator introduces auxiliary variables in order to turn p into a multiaffine polynomial over a larger set of variables. (We elaborate on polarization in Section 3.1.) Formally, following Brändén and Huh [17], for $\kappa \in {\mathbb{Z}}_{\ge 0}^n$, let

and ${\mathbb{R}}_{\kappa }^a[x_{ij}]=\{\text{multi-affine} \text{polynomials} \text{in} \mathbb{R}{[x_{ij}]}_{1\le i\le n,1\le j\le {\kappa }_i}\}.$ The polarization operator ${\mathrm{\Pi }}_{\kappa }:{\mathbb{R}}_{\kappa }[x_i]\to {\mathbb{R}}_{\kappa }^a[x_{ij}]$ replaces every factor $x^{\alpha }$ by

Proposition 3

(Brändén and Huh [17]). If p is d-homogeneous and SLC over ${\mathbb{R}}_{\kappa }[x_i]$, then ${\mathrm{\Pi }}_{\kappa }(p)$ is d-homogeneous and SLC over ${\mathbb{R}}_{\kappa }^a[x_{ij}]$.

We conclude this section by stating a large class of homogeneous polynomials that are known to be strongly log-concave.

Proposition 4

(Brändén and Huh [17]). For $\kappa \in {\mathbb{Z}}_{\ge 0}^n$ and $w:K\to {\mathbb{R}}_{\ge 0}$ a nonnegative weight function, consider

$$g_{\kappa }(x)={\displaystyle \sum _{\alpha \in K}\frac{w(\alpha )}{\alpha !}}{x}^{\alpha },$$(6)

and assume that $g_{\kappa }$ is d-homogeneous. Let $\nu :{\mathbb{Z}}_{\ge 0}^n\to \mathbb{R}\cup \{-\infty \}$ be defined by $\nu (\alpha )=\log (w(\alpha ))$ for $\alpha \in K$ and $\nu (\alpha )=-\infty$ otherwise. If ν is M-concave, then $g_{\kappa }$ is SLC.

2.5. Markov Chains

Let $\mathcal{M}=(\mathrm{\Omega },P)$ be an ergodic and time-reversible Markov chain with state space Ω, transition matrix P, and stationary distribution π (which is unique because of the ergodicity assumption). Reversibility means that $\pi (x)P(x,y)=\pi (y)P(y,x)$ for all $x,y\in \mathrm{\Omega }$. We write $P^t(x,·)$ for the distribution over Ω at time step t with initial state $x\in \mathrm{\Omega }$. The total variation distance $d_{TV}(\pi ,\sigma )$ of two distributions π and σ over Ω is defined as

where, for a distribution σ over Ω, we write $\sigma (S)={\sum }_{x\in S}\sigma (x)$. We say that two distributions π and σ are ϵ-close if $d_{TV}(\pi ,\sigma )\le \mathrm{ϵ}$. The total variation distance of the distribution $P^t(x,·)$ from π at time t with initial state x is denoted by ${\mathrm{\Delta }}_x(t)$. The mixing time of $\mathcal{M}$ with initial state $x\in \mathrm{\Omega }$ is ${\tau }_x(\mathrm{ϵ})=\min \{t:{\mathrm{\Delta }}_x(t^{\prime })\le \mathrm{ϵ} \text{for} \text{all} t^{\prime }\ge t\}.$ Informally, ${\tau }_x(\mathrm{ϵ})$ is the number of steps until the Markov chain is ϵ-close to its stationary distribution given that it is starting in x. A treaty of some more advanced Markov chain notions is given in Appendix A, including the definition of the modified log-Sobolev constant $\rho =\rho (P)$, which can be used to bound the mixing time of a Markov chain. (The definition of ρ is deferred to Appendix A as we actually do not need it in the main body; we only rely on lower bounds on this constant.) It holds that (see, e.g., Bobkov and Tetali [15])

2.5.1. Markov Chain Decomposition.

Let $\mathrm{\Omega }={\mathrm{\Omega }}_1\cup \dots \cup {\mathrm{\Omega }}_m$ be a disjoint partition of the state space Ω. Following Martin and Randall [46], consider $\overline{\pi }(i)=\pi ({\mathrm{\Omega }}_i)={\sum }_{x\in {\mathrm{\Omega }}_i}\pi (x)$ and let $\overline{P}:[m]\times [m]\to [0,1]$ be defined by

The Markov chain on $[m]$ with transition matrix $\overline{P}$ is called the projection chain on the partition ${\{{\mathrm{\Omega }}_i\}}_{i=1,\dots ,m}$. It is time-reversible with respect to the distribution $\overline{\pi }$ over $[m]$. For $i\in [m]$ the restriction chain on Ω_i has transition matrix $P_i:{\mathrm{\Omega }}_i\times {\mathrm{\Omega }}_i\to [0,1]$ given by

It is time-reversible with respect to the distribution ${\pi }_i(x)=\pi (x)/\overline{\pi }(i)$ for $x\in {\mathrm{\Omega }}_i$. In Appendix A.1, we give a Markov chain decomposition result based on the modified log-Sobolev constant ρ.

2.5.2. Base-Exchange Markov Chain.

Let $\mathcal{N}$ be a matroid, and let π be an SLC probability distribution over the set of bases $\mathcal{B}$ given by $\pi (\alpha )\propto w(\alpha )$ for some nonnegative weight function $w:K\to \mathbb{R}$. Here, $\pi (\alpha )\propto w(\alpha )$ means that $\pi (\alpha )=w(\alpha )/({\sum }_{\alpha \in K}w(\alpha ))$. The base-exchange Markov chain on $\mathcal{B}$ is defined by the following transitions, in which $B\in \mathcal{B}$ is the current state of the Markov chain:

• Select an element $e\in B$ uniformly at random and remove it.

• Pick a base $B^{\prime }\in \mathcal{B}$ with $B^{\prime }\supset B-e$ with probability $\propto w(B^{\prime })$ among all such bases $B^{\prime }$.

It is not hard to see, using the base-exchange property, that this procedure defines an ergodic, time-reversible Markov chain with stationary distribution π. Anari et al. [4] show that this chain is rapidly mixing for any matroid $\mathcal{N}$. In particular, in a recent follow-up work, they give a tight mixing time bound (Anari et al. [6]).

Theorem 3

(Anari et al. [6]). Let $\mathcal{N}$ be matroid of rank r, and let π be an SLC probability distribution over the set of bases $\mathcal{B}$ given by $\pi (\alpha )\propto w(\alpha )$ for some weight function $w:\mathcal{B}\to {\mathbb{R}}_{\ge 0}$. Then, the mixing time of the base-exchange random walk satisfies $\tau (\mathrm{ϵ})\le O(r\log (r/\mathrm{ϵ})).$

We note that the mixing time is independent of the size n of the ground set E of the matroid $\mathcal{N}$ as well as the stationary distribution π. In this work, Theorem 3 is essentially only applied to partition and uniform matroids. Furthermore, Cryan et al. [20] show that the modified log-Sobolev of the base-exchange random walk satisfies $\rho \ge 1/r$, where r is the rank of the matroid.

2.6. Sampling Algorithms

Consider a class of (capacitated) congestion games $\mathrm{\Gamma }=(N,E,{({\mathcal{S}}_i)}_{i\in N},{(c_e)}_{e\in E},{(u_e)}_{e\in E})$ with n players and m resources and cost functions $c_e:{\mathbb{Z}}_{\ge 0}\to \mathbb{Q}$ for $e\in E$. Let $w:\mathcal{S}\to {\mathbb{Q}}_{\ge 0}$ be a weight function. In this work, an algorithm for sampling $s\in \mathcal{S}$ according to a distribution ϵ-close to π with $\pi (s)\propto w(s)$ is said to run in (randomized) polynomial time if the number of arithmetic operations can be upper bounded by a polynomial in $n,\hspace{0.17em}m,\log (1/\mathrm{ϵ}),\hspace{0.17em}{\max }_{e,j}\log (c_e(j))$ and ${\max }_s\log (w(s))$. (We assume the strategy sets can be represented in a compact form.) The generation of a uniform random 0/1 bit is considered to be one arithmetic operation.

Remark 3

(Real Numbers). In this work, we use Markov chains whose transition probabilities are, in general, not rational numbers (in particular for the Gibbs distribution). Whenever we use real numbers, it is implicitly assumed that we use sufficiently accurate approximations to these numbers. All our results remain valid when (real-valued) transition probabilities are replaced by sufficiently accurate rational approximations. We note that, roughly speaking, whenever we want to generate Markov chain transitions with probabilities proportional to $e^{-T\mathrm{\Phi }(s)}$ for $s\in \mathcal{S}$, our algorithms run in pseudo-polynomial time in terms of the values of the cost functions of the congestion game under consideration.

All our algorithmic results are based on running Markov chains for a sufficiently long time. We usually write our running time bounds as the product of two factors: the number of steps that we need to run the Markov chain (before it’s close to stationarity) and the complexity of implementing one such step. In particular, in all cases, the transitions probabilities of one step are determined by a sequence of rational numbers $a=(a_1,\dots ,a_z)$ and $q=(q_1,\dots ,q_z)$, and we want to sample an index $i\in [z]$ with probability

Here, z as well as the encoding size of the q_i are $\text{poly}(n,m)$. We refer to $C=C(n,m,a)$ as the computational complexity of sampling an index i according to (approximations of) the preceding probabilities in order not to overload our theorem statements. We say that probabilities of the form (8) are suitable.

2.7. Bipartite Graphs

An (undirected) bipartite graph $G=(A\cup B,F)$ is given by two disjoint sets of nodes $A=\{a_1,\dots ,a_n\}$ and $B=\{b_1,\dots ,b_m\}$ with $F\subseteq \{\{a,b\}:a\in A,b\in B\}$. We say that G has degree sequence $(\mathbf{x},\mathbf{y})$ if $d(a_i)=x_i$ for $i=1,\dots ,n$ and $d(b_i)=y_j$ for $j=1,\dots ,m$, where d(v) denotes the degree of node v in G. We write $\mathcal{G}(\mathbf{x},\mathbf{y})$ for the set of all bipartite graphs on $A\cup B$ with degree sequence $(\mathbf{x},\mathbf{y})$.

3. General Approach

In Section 2.1, we give two possible definitions for the load profile of a strategy profile $s=(s_1,\dots ,s_n)\in {\times }_i{\mathcal{S}}_i$. For general congestion games, we define the resource load profile $\ell (s)=({\ell }_e(s))$ that keeps track of how many players use a particular resource e in s. For symmetric congestion games, we may in addition consider the strategy load profile $z(s)={(z_t(s))}_{t\in {\mathcal{S}}_0}$ that keeps track of how many players use a particular strategy t from the common strategy set $\mathcal{P}$.

A general approach for sampling a strategy profile according to the Gibbs distribution in Sections 4.1 and 6.1 is to first sample a (resource or strategy) load profile α according to approximately the right probability and then sample a strategy profile $s\in \mathcal{S}(\alpha )$ uniformly at random. Remember that the set $\mathcal{S}(\alpha )$ is used to denote all strategy profiles with the given (resource or strategy) load profile. The approach is summarized in Algorithm 1, in which we give the formulation for resource load profiles because we use Rosenthal’s potential $\tilde{\mathrm{\Phi }}:{\mathbb{Z}}_{\ge 0}^m\to \mathbb{R}$ as given in (4). The formulation for strategy load profiles is exactly the same with $\tilde{\mathrm{\Phi }}$ replaced by $\mathrm{\Phi }$ as given in (3). We note that, in order to sample a strategy load profile $s\in \mathcal{S}(\alpha )$ uniformly at random in symmetric games, it suffices to generate a random permutation of the players in $N=\{1,\dots ,n\}$.

Sampling a load profile with the correct probability in our applications corresponds to sampling a base of a discrete polymatroid according to a strongly log-concave distribution. In order to do this, we present a reduction of this problem to that of sampling a base of a matroid according to a strongly log-concave distribution in Section 3.1 (after which we can rely on Theorem 3).

Algorithm 1

(Gibbs Sampler for Congestion Game Γ)

 Input: Congestion game Γ, rationality level $T\ge 0$, and $\mathrm{ϵ}\ge 0$.

 Output: Strategy profile $s\in \mathcal{S}$ according to distribution $\overline{\pi }$ that is ϵ-close to Gibbs distribution π at rationality level T.

• Step I: Sample load profile α according to a distribution ${\sigma }^{\prime }$ that is ϵ-close to ${\pi }^{\prime }$ given by

  $${\pi }^{\prime }(\alpha )=|\mathcal{S}(\alpha )|e^{-T\tilde{\mathrm{\Phi }}(\alpha )}.$$

• Step II: Sample strategy profile $s\in \mathcal{S}(\alpha )$ (approximate) uniformly at random.

3.1. Sampling Bases of Discrete Polymatroids

In this section, we describe how to generate a discrete polymatroid base, according to a strongly log-concave distribution over the set of all polymatroid bases by reducing it to the problem of generating a base of a matroid. This follows more or less directly from Theorem 3 by using the notion of polarization. Polarization can be seen as a functional version of the classic reduction from discrete polymatroids to matroids as given by Helgason [35]. (See also Schrijver [57, chapter 44.6b] for this reduction.)

For a polymatroid $R\subset {\mathbb{Z}}_{\ge 0}^n$, consider a d-homogeneous strongly log-concave polynomial

This is with positive coefficients and support the set of bases ${\mathcal{B}}_R$. Consider the matroid ${\mathcal{N}}_R=(E,\mathcal{I})$ on ground set $E=\{(i,j):1\le i\le n,\hspace{0.17em}1\le j\le d\}$, where $I\in \mathcal{I}$ if and only if the vector $\alpha (I)\in {\mathbb{Z}}_{\ge 0}^n$ given by ${\alpha }_e=|\{f:(e,f)\in I\}|$ satisfies $\alpha (I)\in R$. The fact that ${\mathcal{N}}_R$ is indeed a matroid follows directly from the fact that R is a polymatroid. Note that, for a given $\alpha \in R$, we have

We slightly abuse notation here and write $d=(d,\dots ,d)\in {\mathbb{Z}}^n$ for the all d-vector.

Then, with $\mathrm{\Pi }({\mathcal{B}}_R)$ the set of bases of ${\mathcal{N}}_R$, the polarization $\mathrm{\Pi }(g)$ of g can be written as

where, for $B\in \mathrm{\Pi }({\mathcal{B}}_R)$, we define

Polarization should be interpreted as spreading out the weight $w(\alpha )$ for $\alpha \in R$ equally over all bases $B\in \{A:\alpha (A)=\alpha \}\subseteq \mathrm{\Pi }({\mathcal{B}}_R)$. Proposition 3 implies that $\mathrm{\Pi }(g)$ is also strongly log-concave.

Example 2.

Let $p(x_1,x_2)=x_1^2{x}_2+x_1{x}_2$ so that $\text{supp}(p)=\{(2,0),(1,1)\}$ and take $d=(2,2)$. Then,

$$\begin{array}{ll}{\mathrm{\Pi }}_d(p)\hfill & =\frac{1}{2}{x}_{11}{x}_{12}(x_{21}+x_{22})+\frac{1}{4}(x_{11}+x_{12})(x_{21}+x_{22})\hfill \\ \hfill & =\frac{1}{2}{x}_{11}{x}_{12}{x}_{21}+\frac{1}{2}{x}_{11}{x}_{12}{x}_{22}+\frac{1}{4}{x}_{11}{x}_{21}+\frac{1}{4}{x}_{11}{x}_{22}+\frac{1}{4}{x}_{12}{x}_{21}+\frac{1}{4}{x}_{12}{x}_{22}.\hfill \end{array}$$

Note that, looking at the support of p, we have $\left(\begin{array}{c}d\\ \alpha \end{array}\right)=2$ monomials corresponding to $\alpha =(2,0)$ and $\left(\begin{array}{c}d\\ \alpha \end{array}\right)=4$ monomials corresponding to $\alpha =(1,1)$.

Corollary 1 now follows directly from Theorem 3. It says the following. Suppose the current state of the base-exchange Markov chain after t steps, starting from any state $B_0\in \mathrm{\Pi }({\mathcal{B}}_R)$, is the base $B\in \mathrm{\Pi }({\mathcal{B}}_R)$, and suppose we output the polymatroid base $\alpha (B)$. If t is large enough such that we are in state B with probability close to $w_{\mathrm{\Pi }}(\alpha (B))$ for every $B\in \mathrm{\Pi }({\mathcal{B}}_R)$, then $\alpha (B)$ is outputted with probability close to $w(\alpha (B))$ with $w_{\mathrm{\Pi }}(\alpha (B))$ and $w(\alpha (B))$ as in (10).

Corollary 1.

Let π be the distribution over ${\mathcal{B}}_R$ with $\pi (\alpha )\propto w(\alpha )$, and let ${\mathrm{\Pi }}_{\pi }$ be the distribution over $\mathrm{\Pi }({\mathcal{B}}_R)$ with ${\mathrm{\Pi }}_{\pi }(B)\propto w_{\mathrm{\Pi }}(\alpha (B))$. Let $B\in \mathrm{\Pi }({\mathcal{B}}_R)$, and let ${\mathrm{\Pi }}_{\sigma }^t=P^t(B,·)$ be the distribution over $\mathrm{\Pi }({\mathcal{B}}_R)$ after t steps of the base-exchange Markov chain $\mathcal{M}=(\mathrm{\Pi }({\mathcal{B}}_R),P)$. Let σ^t be the induced distribution over ${\mathcal{B}}_R$ given by ${\sigma }^t(\alpha )={\sum }_{B:\alpha (B)=\alpha }{\mathrm{\Pi }}_{\sigma }^t(B).$ If $d_{TV}({\mathrm{\Pi }}_{\sigma }^t,{\mathrm{\Pi }}_{\pi })\le \mathrm{ϵ}$, then $d_{TV}({\sigma }^t,\pi )\le \mathrm{ϵ}$.

Proof.

We have

$$\begin{array}{ll}2d_{TV}(\sigma ,\pi )\hfill & ={\displaystyle \sum _{\alpha \in {\mathcal{B}}_R}\left|{\displaystyle \sum _{B:\alpha (B)=\alpha }{\sigma }^{\prime }}(B)-w(\alpha )\right|}\hfill \\ \hfill & ={\displaystyle \sum _{\alpha \in {\mathcal{B}}_R}\left|{\displaystyle \sum _{B:\alpha (B)=\alpha }\left[{\sigma }^{\prime }(B)-{\left(\begin{array}{c}d\\ \alpha (B)\end{array}\right)}^{-1}w(\alpha )\right]}\right|}\hspace{0.17em}\hspace{0.17em}(\text{using} (9))\hfill \\ \hfill & \le {\displaystyle \sum _{B\in {\mathcal{B}}^{\prime }}|}{\sigma }^{\prime }(B)-w^{\prime }(\alpha (B))|\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}(\text{triangle} \text{inequality})\hfill \\ \hfill & =2d_{TV}({\sigma }^{\prime },{\pi }^{\prime })\le 2\mathrm{ϵ}.\hfill \end{array}$$

This gives the desired result. □

Remark 4.

It is possible to define a more direct Markov chain on the set of all bases of a given discrete polymatroid and prove that this chain is rapidly mixing (also based on Theorem 3), but this is not needed for our results.

4. Extension Parallel Congestion Games

An extension parallel congestion game is a symmetric congestion game in which the common strategy set $\mathcal{P}$ of the players consists of the o, d-paths in a (directed) extension parallel network $G=(V,A)$ with source o and target d. For two given networks $G_i=(V_i,A_i)$ with source o_i and target d_i for i = 1, 2, let $G^{\prime }=(V_1\cup V_2,A_1\cup A_2)$ be the union of G₁ and G₂. The parallel composition of G₁ and G₂ is the network obtained by identifying o₁ with o₂ and d₁ with d₂. These nodes are the source and target of $G^{\prime }$, respectively. The series composition of G₁ and G₂ is obtained by identifying d₁ with o₂. The node o₁ becomes the source of $G^{\prime }$ and d₂ its target. An extension parallel network consists of (i) a single arc (o, d), (ii) two extension parallel networks in parallel, or (iii) a single arc in series with an extension parallel network. An example is given in Figure 2. For a given extension parallel graph G, we use $\mathcal{P}=\{p_1,\dots ,p_q\}$ to denote all o, d-paths in G. Note that, for an extension parallel network, we have $q\le |A|=m$.

In this section we are always working with strategy load profiles (and so these are sometimes simply referred to as load profiles). The set of all possible strategy load profiles is denoted by $\mathcal{L}=\{\alpha \in {\{1,\dots ,q\}}^n:|\alpha |=n\}$. We consider the potential $\mathrm{\Phi }:\mathcal{L}\to \mathbb{Q}$ defined by $\mathrm{\Phi }(\alpha )=\mathrm{\Phi }(s)$ for some $s\in \mathcal{S}(\alpha )$. This is well-defined as the potential value is the same for any choice of $s\in \mathcal{S}(\alpha )$. The main result that we need in this section is the M-convexity of Rosenthal’s potential for EP congestion games.

Proposition 5

(Fujishige et al. [30]). Let Γ be an extension parallel congestion game. Then, the potential $\mathrm{\Phi }:\mathcal{L}\to \mathbb{Q}$ defined by $\mathrm{\Phi }(\alpha )=\mathrm{\Phi }(s)$ for $s\in \mathcal{S}(\alpha )$ is M-convex.

The M-convexity of $\mathrm{\Phi }$ follows, in a nutshell, from the fact that the collection of sets $\{Q_a:a\in A\}$, where $Q_a=\{p\in \mathcal{L}:a\in p\}$ is the set of all paths containing arc a, form a laminar family in the case of extension parallel networks. This implies that the potential $\mathrm{\Phi }$ is laminar convex, which yields the M-convexity property of the potential. We refer the reader to Fujishige et al. [30] for more details.

Before giving our main result as sketched in Section 1.1, we first give a more direct approach for sampling from the Gibbs distribution in EP congestion games.

4.1. Sampling from the Gibbs Distribution

As mentioned in Section 3 and sketched in Algorithm 1, the high-level algorithmic idea for sampling a strategy profile according to the Gibbs distribution consists of first sampling a load profile α with the correct probability and then a strategy profile from $\mathcal{S}(\alpha )$ uniformly at random. The main result of this section based on this approach is stated in Theorem 4.

Theorem 4.

Let $\mathrm{ϵ}>0$ and $T\ge 0$, and let Γ be an extension parallel congestion game with n players. There is a randomized algorithm $\mathcal{A}$ with output distribution $\overline{\pi }$ over ${\mathcal{P}}^n$ that is ϵ-close to the Gibbs distribution π at rationality level T and runs in (expected) time $O(C·n\log (n/\mathrm{ϵ}))$ with C the complexity of implementing one step of a base-exchange Markov chain with suitable probabilities (see Section 2.6).

Proof.

Note that, for an extension parallel congestion game, the number of strategy profiles corresponding to a given load profile α is $|\mathcal{S}(\alpha )|=n!/\alpha !$. (This is the number of ways in which we can assign n labeled balls to bins $b_1,\dots ,b_q$, where b_i contains α_i balls.)

Lemma 1.

The n-homogeneous generating polynomial

$$g(x_1,\dots ,x_n)={\displaystyle \sum _{\alpha \in {[q]}^n:|\alpha |=n}|}\mathcal{S}(\alpha )|e^{-T\mathrm{\Phi }(\alpha )}{x}^{\alpha }={\displaystyle \sum _{\alpha \in {[q]}^n:|\alpha |=n}\frac{n!}{\alpha !}}{e}^{-T\mathrm{\Phi }(\alpha )}{x}^{\alpha }$$(11)

is strongly log-concave. Hence, the distribution ${\pi }^{\prime }$ over $\mathcal{L}$ given by ${\pi }^{\prime }(\alpha )\propto \frac{n!}{\alpha !}{e}^{-T\mathrm{\Phi }(\alpha )}$ for $\alpha \in \mathcal{L}$ is strongly log-concave.

Proof.

Strong log-concavity is preserved under scalar multiplication by Proposition 2, so it suffices to show that

$$\frac{1}{n!}g(x)={\displaystyle \sum _{\alpha \in {[q]}^n:|\alpha |=n}\frac{e^{-T\mathrm{\Phi }(\alpha )}}{\alpha !}}{x}^{\alpha }$$

is strongly log-concave. In turn, by Proposition 4, it is sufficient to show that $\log \left(e^{-T\mathrm{\Phi }(\alpha )}\right)=-T\mathrm{\Phi }(\alpha )$ is an M-concave function on its effective domain. As $T\ge 0$, this is equivalent to showing that $\mathrm{\Phi }(\alpha )$ is M-convex on its effective domain $\mathcal{L}=\{\alpha \in {[q]}^n:|\alpha |=n\}$. This follows from Proposition 5. □

Because of Lemma 1, the polarization $\mathrm{\Pi }(g)$ of g in (11) is also strongly log-concave. The support of $\mathrm{\Pi }(g)$ can be seen as the bases of the n-uniform matroid $\mathcal{N}$ on ground set $\{(i,j):1\le i,j\le n\}$. Our algorithm now consists of first running the base-exchange Markov chain for $O(n\log (n/\mathrm{ϵ}))$ steps, starting from any initial base. We output $\alpha (B)$, where B is the state we are in after the $O(n\log (n/\mathrm{ϵ}))$ steps that were carried out. The resulting distribution ${\sigma }^{\prime }(\alpha )$ over $\mathcal{L}$ satisfies $d_{TV}({\sigma }^{\prime },{\pi }^{\prime })\le \mathrm{ϵ}$ by Corollary 1. We then uniformly at random choose a strategy profile from $\mathcal{S}(\alpha )$. Let $\overline{\pi }$ be the resulting output distribution over $\mathcal{S}$. It remains to show that $d_{TV}(\pi ,\overline{\pi })\le \mathrm{ϵ}$, which the following calculation shows. Note that

$$\overline{\pi }(s)=\frac{\alpha !}{n!}{\sigma }^{\prime }(\alpha ),$$

where $\alpha =\ell (s)$ is the load profile corresponding to strategy s. Then,

$$\begin{array}{ll}{\displaystyle \sum _{s\in \mathcal{S}}|}\overline{\pi }(s)-\pi (s)|\hfill & ={\displaystyle \sum _{\alpha }{\displaystyle \sum _{s\in \mathcal{S}:\ell (s)=\alpha }|}}\overline{\pi }(s)-\pi (s)|\hfill \\ \hfill & ={\displaystyle \sum _{\alpha }{\displaystyle \sum _{s\in \mathcal{S}:\ell (s)=\alpha }\left|\frac{\alpha !}{n!}{\sigma }^{\prime }(\alpha )-e^{-T\mathrm{\Phi }(\alpha )}\right|}}\hfill \\ \hfill & ={\displaystyle \sum _{\alpha }\frac{n!}{\alpha !}}\left|\frac{\alpha !}{n!}{\sigma }^{\prime }(\alpha )-e^{-T\mathrm{\Phi }(\alpha )}\right|\hfill \\ \hfill & ={\displaystyle \sum _{\alpha }\left|{\sigma }^{\prime }(\alpha )-\frac{n!}{\alpha !}{e}^{-T\mathrm{\Phi }(\alpha )}\right|}\hfill \\ \hfill & \le 2\mathrm{ϵ}.\hfill \end{array}$$

We conclude with analyzing the running time of the algorithm. One step of the base-exchange Markov chain can be implemented in time O(C) by definition. Generating an $s\in \mathcal{S}(\alpha )$ uniformly at random can be done by generating a uniform random permutation μ of $\{1,\dots ,n\}$. We set $s_i=p_1$ for players $i=\mu (1),\dots ,\mu ({\alpha }_1),\hspace{0.17em}{s}_i=p_2$ for players $\mu ({\alpha }_1+1),\dots ,\mu ({\alpha }_2+1)$, and so on. Generating a uniform random permutation can be done in time $O(n\log (n))$. □

4.2. Relaxed Logit Dynamics

In this section, we give the proof of Theorem 1, reformulated in Theorem 5. Formally, the relaxed logit dynamics Markov chain with current state $s\in {\mathcal{P}}^n$ proceeds by

• With probability 1/2: Select two players $i,j\in N$ uniformly at random and transition to $s^{\prime }$ given by

  $$s_k^{\prime }=\{\begin{array}{ll}{s}_i\hfill & \text{if} k=j\hfill \\ s_j\hfill & \text{if} k=i\hfill \\ s_k\hfill & \text{otherwise}.\hfill \end{array}$$

• With probability 1/2: Perform a transition according to the logit dynamics (as in Section 2.2).

We note that, for any symmetric congestion game, this is a well-defined ergodic, time-reversible Markov chain with the Gibbs distribution as stationary distribution.

Theorem 5.

For an extension parallel congestion game Γ with common strategy set $\mathcal{P}$ and initial state $s\in {\mathcal{P}}^n$, the mixing time of the relaxed logit dynamics Markov chain at rationality level $T\ge 0$ satisfies 

$${\tau }_s(\mathrm{ϵ})\le n^3\left(\log n+\log \log |\mathcal{P}|+\log \left(\frac{2T{\mathrm{\Phi }}_{\max }}{{\mathrm{ϵ}}^2}\right)\right),$$

where ${\mathrm{\Phi }}_{\mathit{max}}={\max }_{r\in \mathcal{S}}\mathrm{\Phi }(r)$ is the maximum value attained by Rosenthal’s potential over ${\mathcal{P}}^n$.

Compared with the mixing time of the (nonrelaxed) logit dynamics for general games (Auletta et al. [12]), we get a doubly exponential improvement in terms of the dependence on $T{\mathrm{\Phi }}_{\mathit{max}}$ (at the cost of a small polynomial increase in the dependence on n).

The proof of Theorem 5 relies on a Markov chain decomposition argument. The idea underlying Markov chain decomposition is to divide the state space into a disjoint union of smaller sets. If one can show that the Markov chain is rapidly mixing when restricted to every smaller set and that it is easy to transition between the different smaller sets, then the original chain is rapidly mixing as well. In our setting, every load profile is identified with a smaller set, which is formed by all strategy profiles inducing the given load profile. To prove rapid mixing of the smaller sets, we compare the restricted chains to the so-called random transposition Markov chain, which is known to be rapidly mixing; see, for example, Goel [32]. To prove rapid mixing between the smaller sets, that is, between load profiles, we use again the connection to the base-exchange Markov chain but in a different way than in the proof of Theorem 4. In order to prove rapid mixing of all the Markov chains involved, we show that their modified log-Sobolev constant can be bounded appropriately. Recall that this is a quantity that can be used to upper bound the mixing time of a Markov chain (see Appendix A).

Proof of Theorem 5.

We use a Markov chain decomposition argument based on the two operations that define the relaxed logit dynamics Markov chain. We first partition the state space $\mathcal{S}={\mathcal{P}}^n$ naturally based on load profiles by setting ${\mathrm{\Omega }}_{\alpha }=\mathcal{S}(\alpha )$ for $\alpha \in \mathcal{L}$, whereas before we have $\mathcal{L}=\{\alpha :\alpha \in {[q]}^n \text{and} \text{|}\alpha \text{|}=n\}$. Our proof approach is to apply the Markov chain decomposition theorem of Hermon and Salez [36] as given in Theorem A.1. In particular, for this, we need to bound the modified log-Sobolev constants of the projection and restriction chains. We start with the modified log-Sobolev constant $\overline{\rho }$ of the projection chain.

The projection chain $\overline{P}$ has state space $\mathcal{L}$ and stationary distribution $\overline{\pi }(\alpha )=|\mathcal{S}(\alpha )|e^{-T\mathrm{\Phi }(\alpha )}$ for $\alpha \in \mathcal{L}$. Let $\alpha ,\beta \in \mathcal{L}$ such that ${\sum }_e|{\alpha }_e-{\beta }_e|=2$; that is, there exist paths p and $p^{\prime }$ such that

$${\alpha }_e=\{\begin{array}{ll}{\beta }_e+1\hfill & \text{if} e=p\hfill \\ {\beta }_e-1\hfill & \text{if} e=p^{\prime }\hfill \\ {\alpha }_e\hfill & \text{if} e\in E\backslash \{p,p^{\prime }\}.\hfill \end{array}$$

In this case, we say that α and β are adjacent load profiles differing on paths p and $p^{\prime }$. Note that, if $s\in \mathcal{S}(\alpha )$ and $s^{\prime }\in \mathcal{S}(\beta )$ are such that they differ by a deviation of some player i from path p to path $p^{\prime }$, then $P(s,s^{\prime })=\frac{1}{2n}\frac{\exp (-T\mathrm{\Phi }(p^{\prime },s_{-i}))}{Z^{\prime }}$, and this expression is the same for every such player i with s_i = p. Moreover, note that $\pi (x)=\pi (y)$ for any two strategy profiles $x,y\in \mathcal{S}(\alpha )$.

For some fixed choice of strategy profile $x\in {\mathrm{\Omega }}_{\alpha }$ and player i using path p, that is, s_i = p, the transition probabilities for adjacent load profiles can then be seen to equal (with Z and $Z^{\prime }$ are the normalizing constants as in Section 2.2)

$$\begin{array}{ll}2\overline{P}(\alpha ,\beta )=\frac{1}{\overline{\pi }(\alpha )}{\displaystyle \sum _{x\in {\mathrm{\Omega }}_{\alpha },y\in {\mathrm{\Omega }}_{\beta }}\pi }(x)P(x,y)\hfill & =\frac{\pi (x)}{\overline{\pi }(\alpha )}\frac{{\alpha }_p}{n}|\mathcal{S}(\alpha )|\frac{\exp (-T\mathrm{\Phi }(\beta ))}{Z^{\prime }}\hfill \\ \hfill & =\frac{\exp (-T\mathrm{\Phi }(\alpha ))/Z}{|\mathcal{S}(\alpha )|\exp (-T\mathrm{\Phi }(\alpha ))/Z}\frac{{\alpha }_p}{n}|\mathcal{S}(\alpha )|\frac{\exp (-T\mathrm{\Phi }(\beta ))}{Z^{\prime }}\hfill \\ \hfill & =\frac{{\alpha }_p}{n}\frac{\exp (-T\mathrm{\Phi }(\beta ))}{Z^{\prime }}\hfill \\ \hfill & =2{\alpha }_{p}P(s,s^{\prime }),\hfill \end{array}$$(12)

where the last equality is true for any choice of $s\in \mathcal{S}(\alpha )$ and $s^{\prime }\in \mathcal{S}(\beta )$. Note that this implies that, for any $\alpha ,\beta \in \mathcal{L},\hspace{0.17em}s\in \mathcal{S}(\alpha )$ and $s^{\prime }\in \mathcal{S}(\beta )$, we have

$$\frac{P(s,s^{\prime })}{\overline{P}(\alpha ,\beta )}=\frac{1}{{\alpha }_p}\ge \frac{1}{n}.$$(13)

The lower bound of $1/n$ serves as our lower bound on χ as defined in Appendix A.1. In order to bound the modified log-Sobolev constant of the projection chain, one can use a comparison argument (as defined in Appendix A.2) with the base-exchange Markov chain on the support of the polarization $\mathrm{\Pi }(g_{\mathrm{\Gamma }})$ of $g_{\mathrm{\Gamma }}$ as in (11). In this section, the support corresponds to the set of bases of an n-uniform matroid. In particular, it holds that

$$\overline{\rho }\ge \frac{1}{n}·\rho (\mathrm{\Pi }(\mathcal{L}))\ge \frac{1}{n^2},$$(14)

where $\rho (\mathrm{\Pi }(\mathcal{L}))$ is the modified log-Sobolev constant of the base-exchange Markov chain on the support of $\mathrm{\Pi }(g)$ with g as in (11). The second inequality comes from the fact that $\rho (\mathrm{\Pi }(\mathcal{L}))\ge 1/n$ as is shown by Cryan et al. [20]. The first inequality is somewhat tedious to prove as it requires a Markov chain comparison argument between two Markov chains on different state spaces and is, therefore, deferred to Appendix B.

We continue with bounding the modified log-Sobolev constant of the restriction chains. In order to do this, we use a comparison argument with the random transposition Markov chain on the set S_k of all permutations of $\{1,\dots ,k\}$. Given a permutation σ, this chain proceeds by selecting two positions a and b uniformly at random and interchanging the positions of the elements $\sigma (a)$ and $\sigma (b)$. With ρ_rt denoting the modified log-Sobolev constant of this chain, it follows that, for every $\alpha \in \mathcal{L}$, we have

$${\rho }_{\alpha }\ge {\rho }_{rt}\ge \frac{1}{n-1}$$(15)

using the fact that ${\rho }_{rt}\ge 1/(n-1)$ as shown by Goel [32]. This comparison argument is also deferred to Appendix B. Now, applying Theorem A.2, it follows that $\overline{\rho }\ge 1/n^3$. Plugging this into (7), it then follows that

$${\tau }_s(\mathrm{ϵ})\le n^3\left(\log n+\log \log |\mathcal{P}|+\log \left(\frac{2T{\mathrm{\Phi }}_{\max }}{{\mathrm{ϵ}}^2}\right)\right)$$

using that $\pi {(s)}^{-1}\le |\mathcal{P}{|}^n{e}^{-T{\mathrm{\Phi }}_{\max }}$ for every $s\in \mathcal{S}$ because of the nonnegativity of the cost functions. □

4.3. Uniform Sampling of Pure Nash Equilibria

In Theorem 6, we show that the result in Theorem 4 also implies that, for an extension parallel congestion game Γ, we can (approximate) uniformly at random sample a pure Nash equilibrium from the set $\text{NE}(\mathrm{\Gamma })$ of all pure Nash equilibria of Γ in pseudo-polynomial time. That is, we sample every $s\in \text{NE}(\mathrm{\Gamma })$ with probability $\approx 1/|\text{NE}(\mathrm{\Gamma })|$. The (approximate) uniform sampling of combinatorial objects has received a lot of attention in the last 30 years, in particular within the area of theoretical computer science. However, to the best of our knowledge, no nontrivial results for (pure) Nash equilibria are known despite the fact that the problem of computing Nash equilibria has received much attention.

For the proof of Theorem 6, we use the fact that Nash equilibria are precisely the strategy profiles minimizing Rosenthal’s potential in EP congestion games. Furthermore, we exploit the fact that Theorem 4 holds for any rationality level $T\ge 0$. In particular, if we set T large enough, then most weight in the stationary distribution is assigned to profiles minimizing Rosenthal’s potential (under the assumption that the cost functions are integer-valued). This means that, with high probability, Algorithm 1 outputs a strategy profile minimizing Rosenthal’s potential with domain $\mathcal{S}$. (Whenever we refer to Algorithm 1, we mean the implementation of the high-level approach as given in the previous section.) We use the Gibbs distribution with base two instead of e to avoid having to work with real numbers.

Theorem 6.

Let $\mathrm{ϵ}>0$, and let Γ be an extension parallel congestion game with integer-valued cost functions and n players. There is a randomized algorithm $\mathcal{A}$ with output distribution $\overline{\pi }$ over $\text{NE}(\mathrm{\Gamma })$ that is ϵ-close to the uniform distribution over $\text{NE}(\mathrm{\Gamma })$ and runs in (expected) time polynomial in $n,m,{\mathrm{\Phi }}_{\max }$ and $\log (1/\mathrm{ϵ})$, where ${\mathrm{\Phi }}_{\max }$ is the maximum value attained by Rosenthal’s potential.

For the proof of Theorem 6, we use the following correspondence between Nash equilibria and strategy profiles minimizing Rosenthal’s potential.

Proposition 6

(Fotakis [29], Holzman and Law-Yone [37]). The set of strategy profiles $\text{NE}(\mathrm{\Gamma })$ of an extension parallel congestion game Γ coincides with the set of strategy profiles that minimize Rosenthal’s potential as in (3).

Proof of Theorem 6.

We first show that, for T sufficiently large in the algorithm used to prove Theorem 4, most weight is assigned to strategy profiles minimizing Rosenthal’s potential. We apply the idea in Algorithm 1 used to prove Theorem 4 with base two instead of base e. Remember that q is the number of (o, d)-paths in the extension parallel network of the game Γ. Let $\varphi =\mathrm{\Phi }(s)$ be the common potential value of all strategy profiles $s\in \text{NE}(\mathrm{\Gamma })$. For any other strategy profile $s^{\prime }\in \mathcal{S}\backslash \text{NE}(\mathrm{\Gamma })$, we have

$$2^{-T\mathrm{\Phi }(s^{\prime })}\le 2^{-T(\varphi +1)}=2^{-T}{e}^{-T\varphi }$$

by assumption that all cost functions are integer-valued. As there are q strategies to choose from for every player, we have $|\mathcal{S}|=q^n=2^{n{\log }_2(q)}$. This implies that the Gibbs distribution π over $\mathcal{S}$ with rationality level $T=\lceil n{\log }_2(q)+{\log }_2(2/\mathrm{ϵ})\rceil$ satisfies

$$\pi (\mathcal{S}\backslash \text{NE}(\mathrm{\Gamma }))={\displaystyle \sum _{s\in \text{NE}(\mathrm{\Gamma })}{2}^{-T\mathrm{\Phi }(s^{\prime })}}\le 2^{n\log (q)}{2}^{-T}{2}^{-T\varphi }\le \frac{\mathrm{ϵ}}{2}·\pi (\text{NE}(\mathrm{\Gamma })).$$(16)

The algorithm for sampling an (almost) uniform sample from $\text{NE}(\mathrm{\Gamma })$ now works as follows. First, compute a strategy profile minimizing Rosenthal’s potential in order to determine $\varphi$. This can be done efficiently; see, for example, Fotakis [29]. Then, run Algorithm 1 with $T=\lceil n{\log }_2(q)+{\log }_2(2/\mathrm{ϵ})\rceil$ and ${\mathrm{ϵ}}^{\prime }=\mathrm{ϵ}/2$. If the resulting strategy profile has potential value $\varphi$, output this strategy profile, and otherwise, rerun Algorithm 1 until it does. Note that, with probability at least $(1-\mathrm{ϵ}/2)$, Algorithm 1 outputs a strategy profile with potential value $\varphi$ in one run. A simple argument then shows that the output distribution is ϵ-close to the uniform distribution over $\text{NE}(\mathrm{\Gamma })$ as desired. □

Remark 5.

The pseudo-polynomial dependence coming from the polynomial dependence on ${\mathrm{\Phi }}_{\max }$ rather than ${\log }_2({\mathrm{\Phi }}_{\max })$ arises from the fact that we have to compute transition probabilities of the form $2^{a_i}/{\sum }_i{2}^{a_i}$, where the a_i are integers, which requires $\mathrm{\Omega }\left({\sum }_i{a}_i\right)$ random 0/1 bits (following the notion of suitable probabilities in Section 2.6). However, there is no a priori reason that the problem of (approximately) sampling pure Nash equilibria according to the uniform distribution requires pseudo-polynomial (in the input size of the original congestion game) time as opposed to sampling from the Gibbs distribution. We leave open the question of finding a (truly) polynomial time algorithm.

5. Max-k-Cut Games

In this section, we give an application of Theorem 5 to max-k-cut games. We first repeat their definition and show how they can be modeled as a special case of extension-parallel congestion games in case the graph G of the game is complete.

A max-k-cut game is given by an undirected graph $G=(V,E)$ whose nodes V are the players. Every $v\in V$ has strategy set $[k]=\{1,\dots ,k\}$ whose elements are referred to as colors. With $\mathcal{N}(i)=\{j\in V:\{i,j\}\in E\}$ the set of neighbors of node $i\in V$, the utility of player i is the number of neighbors $U_i(s)=\{j\in \mathcal{N}(i):s_i\ne s_j\}$ that choose a different color (hence, the name max-k-cut game, referring to the well-known max-cut problem (Garey and Johnson [31])). Equivalently, the cost of a player is given by the number of neighbors that choose the same color; that is, for $s\in {[k]}^n$, we have $C_i(s)=|\{j\in \mathcal{N}(i):s_i=s_j\}|,$ and in terms of (relaxed) logit dynamics, considering the cost of a player is equivalent to considering its utility because $C_i(s)+U_i(s)=|\mathcal{N}(i)|$ is independent of the chosen color in s.

When $G=(V,E)$ is the complete graph, a max-k-cut game can naturally be modeled as a so-called symmetric singleton congestion game in which $E=[k]$ is the set of resources and every resource is equipped with the cost function $c_e(x)=x-1$ for $e\in [k]$. Indeed, if $\ell$ players are using resource/color e, then the cost for every player choosing color $\ell$ is $\ell -1$, which is precisely the number of neighbors choosing that color for every such node (because of the fact that G is the complete graph). Such games are easily seen to be special cases of extension parallel congestion games.

The relaxed logit dynamics may then be interpreted as the process by which either a player changes its color (according to the transitions as specified by the logit dynamics) or the colors of two players are interchanged. As ${\mathrm{\Phi }}_{\max }\le k(k-1)/2$, which is the total number of edges of the graph G, we obtain the following corollary of Theorem 5.

Corollary 2.

For a max-k-cut game played on the complete graph $G=(V,E)$ with initial state $s\in {[k]}^n$, where $n=|V|$, the mixing time of the relaxed logit dynamics Markov chain at rationality level $T\ge 0$ satisfies

$${\tau }_s(\mathrm{ϵ})\le n^3\left(\log n+\log \log k+\log \left(\frac{Tk(k-1)}{{\mathrm{ϵ}}^2}\right)\right).$$

Levin et al. [44] show that the mixing time of the (nonrelaxed) logit dynamics on the complete graph depends on the parameter $T\ge 0$ in the case k = 2. They show there is a critical value T_c so that, when $T\le T_c$, the logit dynamics converge quickly to the Gibbs distribution, whereas when $T>T_c$, the dynamics converge slowly. What Corollary 2 shows is that, if one in addition is allowed to randomly interchange the colors of two nodes, the slow mixing for $T>T_c$ can be resolved. We believe this to be of independent interest.

6. Capacitated Uniform Congestion Games

In this section, we consider u-capacitated k-uniform congestion games for vectors $u=(u_1,\dots ,u_m)$ and $k=(k_1,\dots ,k_n)$. We write $K=|k|={\sum }_i{k}_i$ and $U=|u|={\sum }_e{u}_e$. The vector u models the capacities of the resources $e\in E$, that is, the variables ${(u_e)}_{e\in E}$ as defined in Section 2.1. The strategy set of player $i\in N$ is given by all subsets $S\subseteq E$ of cardinality $|S|=k_i$, that is, the bases of the k_i-uniform matroid on E. We write $\mathrm{\Gamma }(u,k)$ for the collection of all u-capacitated k-uniform congestion games. We remark that, in this section, load profiles refer to resource load profiles as defined in Section 2.1 and no longer to path load profiles as considered in Section 4.

Note that we can naturally model a feasible strategy profile in $s=(s_1,\dots ,s_n)\in \mathcal{S}$ of a capacitated uniform congestion game as a (simple) bipartite graph $G=(N\cup E,F)\in \mathcal{G}(u,k)$ on $N\cup E$: there is an edge $\{i,e\}\in F$ if and only if player $i\in N$ uses resource $e\in E$ in its strategy s_i. The main result needed in this section is stated in Proposition 7. We use the notation ${[x]}_b=x(x-1)\cdots (x-b+1)$ for $x,b\in {\mathbb{Z}}_{\ge 1}$. For a bipartite degree sequence $(k,\alpha )$, we then write $K_b={\sum }_{i=1}^n\hspace{0.17em}{[k_i]}_b$ and $A_b={\sum }_{j=1}^m\hspace{0.17em}{[{\alpha }_j]}_b$. Note that $K=A=K_1=A_1$.

Proposition 7

(McKay [49]). Let $\mathcal{D}$ be the collection of all bipartite degree sequences $(k,\alpha )$ for which $1\le k_{\max }{\alpha }_{\max }=o\left(K^{1/4}\right)$. Then,

$$|\mathcal{G}(k,\alpha )|=\frac{K!}{{\prod }_i{k}_i!{\prod }_j{\alpha }_j!}\exp \left(-\frac{K_2}{K^2}·A_2+O\left(\max {\{k_{\max },{\alpha }_{\max }\}}^4/K\right)\right)$$

as $K\to \infty$.

6.1. Sampling from the Gibbs Distribution

In this section, we give an (almost) polynomial time sampling algorithm that samples from a distribution that is close to the Gibbs distribution provided the game is sufficiently large. That is, we show that, for a large class of pairs (u, k), we can sample from a distribution close to the Gibbs distribution.

We follow again the high-level approach in Algorithm 1. The set of all feasible load profiles is now given by $L(k,u)=\left\{\alpha :0\le \alpha \le u \text{and} |\alpha |={\sum }_i{k}_i\right\}$ and $\mathcal{S}(\alpha )=\mathcal{G}(k,\alpha )$ for any feasible load profile $\alpha \in L(k,u)$. Recall that we want to sample an $\alpha \in L(k,u)$ with probability proportional to (approximately) $\approx |\mathcal{S}(\alpha )|e^{-T\mathrm{\Phi }(\alpha )},$ and then sample a strategy profile $s\in \mathcal{S}(\alpha )$ with probability $\approx 1/|\mathcal{S}(\alpha )|$.

A couple of problems arise here compared with the case of extension parallel congestion games. First of all, there is no polynomial time algorithm known to compute the numbers $w_{\alpha }=|\mathcal{S}(\alpha )|$. (In fact, it is still an open question whether this problem is $\#$ P-complete; see, e.g., Jerrum et al. [40].) Instead, we rely on a fully polynomial randomized approximation scheme for computing approximations ${\widehat{w}}_{\alpha }$ to the numbers $w_{\alpha }$ up to arbitrary precision (Bezáková et al. [13], Jerrum et al. [40]). Second, in this case, the polynomial

is in general not strongly log-concave. We overcome this problem by showing that we can restore strong log-concavity “approximately” when the game becomes large and when there are, in addition, suitable capacity constraints. We do this by using asymptotic enumeration formulas for the number of bipartite graphs with a given degree sequence, an area that has received considerable attention in combinatorics. It turns out that replacing $|\mathcal{S}(\alpha )|$ by an asymptotic approximation $\varphi (\alpha )$ in (17) gives rise to a strongly log-concave polynomial.

Finally, the problem of sampling a strategy profile $s\in \mathcal{S}(\alpha )=\mathcal{G}(k,\alpha )$ now corresponds to that of sampling a bipartite graph with degree sequence $(k,\alpha )$ for which many algorithms are known. The main result of this section is given in Theorem 7.

Theorem 7.

Let $\mathrm{ϵ}\ge 0$, let π be the Gibbs distribution at rationality level $T\ge 0$, and let $\mathcal{D}$ be the class of all congestion games $\mathrm{\Gamma }(k,u)$ satisfying

$$1\le k_{\max }{u}_{\max }=o\left(K^{1/4}\right).$$(18)

There is a randomized algorithm $\mathcal{A}$ for the class $\mathcal{D}$ and a constant $K_0\ge 0$ such that the output distribution $\overline{\sigma }$ over $\mathcal{S}$ has the property that $d_{TV}(\overline{\sigma },\pi )\le \mathrm{ϵ}$ whenever $K\ge K_0$. The algorithm runs in (expected) time

$$C·n\left(\log n+\log \log |\mathcal{P}|+\log \left(\frac{2T{\mathrm{\Phi }}_{\max }}{{\mathrm{ϵ}}^2}\right)\right)$$

with $C(n,m,\mathrm{ϵ},{\mathrm{\Phi }}_{\max })= \text{poly}(1/\mathrm{ϵ},n,m,{\mathrm{\Phi }}_{\max })$.

Proof.

Setting

$$\varphi (\alpha )=\frac{K!}{k!\alpha !}\exp \left(-\frac{K_2}{K^2}·A_2\right),$$

it follows, assuming (18) holds, that, for any $0\le \alpha \le u$, we have $\varphi (\alpha )=(1+o(1))|\mathcal{S}(\alpha )|$, where o(1) is with respect to $K\to \infty$. In particular, if $K\ge K_0$ for K₀ large enough, it follows that

$$\frac{1}{2}|\mathcal{S}(\alpha )|\le \varphi (\alpha )\le \frac{3}{2}|\mathcal{S}(\alpha )|.$$(19)

The next step is now to show that replacing $|\mathcal{S}(\alpha )|$ by $\varphi (\alpha )$ in (20) gives rise to a strongly log-concave polynomial. The crucial observation here is to see that $A_2={\sum }_j{\alpha }_j({\alpha }_j-1)$ is a separable convex function.

Lemma 2.

The K-homogeneous generating polynomial

$$g(x_1,\dots ,x_n)={\displaystyle \sum _{\alpha \in L(k,u)}\frac{K!}{k!\alpha !}}·\exp \left(-\frac{K_2}{K^2}·{\displaystyle \sum _{j=1}^m{\alpha }_j}({\alpha }_j-1)\right)\exp \left(-T\mathrm{\Phi }(\alpha )\right)·x^{\alpha }$$(20)

is strongly log-concave.

Proof.

Following the proof of Lemma 1, first observe that (and note that k and all quantities involving k are considered fixed)

$$\frac{k!}{K!}g(x_1,\dots ,x_n)={\displaystyle \sum _{\alpha \in L(k,u)}\frac{1}{\alpha !}}\exp \left(-\frac{K_2}{K^2}·{\displaystyle \sum _{j=1}^m{\alpha }_j}({\alpha }_j-1)\right)\exp \left(-T\mathrm{\Phi }(\alpha )\right)x^{\alpha }$$

is strongly log-concave as well because of Proposition 2. Then, in order to apply Proposition 4, it suffices to show that

$$-\left(\frac{K_2}{K^2}·{\displaystyle \sum _{j=1}^m{\alpha }_j}({\alpha }_j-1)+T\mathrm{\Phi }(\alpha )\right)$$

is M-concave over its domain L(k, u). (Formally speaking, we define it to be $-\infty$ outside of L(k, u).) This follows directly from the fact that both ${\sum }_{j=1}^m{\alpha }_j({\alpha }_j-1)$ and $\mathrm{\Phi }(\alpha )$ are separable convex functions in $({\alpha }_1,\dots ,{\alpha }_m)$ over the (effective) domain L(k, u), and the fact that $K_2,K,T\ge 0$. Separable convex functions (over effective domain L(k, u)) are known to be M-convex (Murota [52]). □

One can now carry out similar steps as in the proof of Theorem 4, albeit with some modifications. Again, the polarization $\mathrm{\Pi }(g)$ of g as in (20) is strongly log-concave as well. The support of $\mathrm{\Pi }(g)$ is now the set of bases of the n-truncation of a partition matroid $\mathcal{N}$ on ground set $E={\cup }_j{E}_j$, where $E_j=\{(j,i):1\le i\le n\}$ with $A\subseteq E$ is independent if and only if $|A\cap E_j|\le u_j$ for $j=1,\dots ,m$. It follows that the modified log-Sobolev constant of this chain satisfies ${\rho }_{\mathcal{N}}\ge 1/n$. Using a Markov chain comparison argument as described in Remark A.1 in Appendix A.2 with $\delta =1/2$ because of (19) then yields that the modified log-Solev constant ρ of the Markov chain in which we use the original quantities $|\mathcal{S}(\alpha )|$, instead of the approximation $\varphi (\alpha )$, satisfies $\rho \ge 1/(3n)$.

The algorithmic problem is now that we cannot compute the quantities $w_{\alpha }=|\mathcal{S}(\alpha )|$ exactly in polynomial time, so one step of this base-exchange Markov chain cannot be implemented efficiently. Nevertheless we can use the approximation scheme of Bezáková et al. [13] to compute approximations ${\widehat{w}}_{\alpha }$ up to arbitrary precision in time polynomial in n, m and $1/\mathrm{ϵ}$ (this gives the dependence of $1/\mathrm{ϵ}$ in C). A Markov chain comparison argument (again as in Remark A.1) then implies that the base-exchange Markov chain on $\mathcal{N}$ using these approximations is also rapidly mixing, and in particular, it is sufficient to run the chain

steps and then output $\alpha (B)$, where B is the current base after having run the chain for the aforementioned number of steps.

We conclude with the sampling of a strategy profile from the set $\mathcal{S}(\alpha )$ uniformly at random, which now requires sampling a bipartite graph from the set $\mathcal{G}(k,\alpha )$ uniformly at random. One algorithm to do this for degree sequences satisfying the condition in Proposition 7 is that of Arman et al. [7] that runs in expected polynomial time. For an overview of algorithms that can be used to (approximately) sample a bipartite graph with a given degree sequence, see, for example, Dyer et al. [23]. □

Remark 6.

We remark here that the algorithm of Theorem 7 does not run in polynomial time as described in Section 2.6 because of the dependence of C on $1/\mathrm{ϵ}$ (as opposed to the required $\log (1/\mathrm{ϵ})$). This dependence arises because of the algorithmic approximations of the numbers $|\mathcal{S}(\alpha )|$ used in the proof of Theorem 7. Alternatively, we could just use the approximations $\varphi (\alpha )$ straight away. However, these predictions only become accurate when $K\to \infty$, so this gives a weaker result in terms of closeness to the Gibbs distribution.