Deep Learning for Computing Convergence Rates of Markov Chains

Authors: yanlin  qu
jose  blanchet
peter  glynn
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 05:35:23

Abstract: Convergence rate analysis for general state-space Markov chains is fundamentally important in operations research (stochastic systems) and machine learning (stochastic optimization). This problem, however, is notoriously difficult because traditional analytical methods often do not generate practically useful convergence bounds for realistic Markov chains. We propose the Deep Contractive Drift Calculator (DCDC), the first general-purpose sample-based algorithm for bounding the convergence of Markov chains to stationarity in Wasserstein distance. The DCDC has two components. First, inspired by the new convergence analysis framework in (Qu et.al, 2023), we introduce the Contractive Drift Equation (CDE), the solution of which leads to an explicit convergence bound. Second, we develop an efficient neural-network-based CDE solver. Equipped with these two components, DCDC solves the CDE and converts the solution into a convergence bound. We analyze the sample complexity of the algorithm and further demonstrate the effectiveness of the DCDC by generating convergence bounds for realistic Markov chains arising from stochastic processing networks as well as constant step-size stochastic optimization.

Main paper: 
 introduction: 
 General state-space Markov chains are indispensable in a wide array of fields due to their flexibility and applicability in modeling random dynamical systems. To analyze the long-term behavior of these Markovian models, estimating the rate of convergence to equilibrium is critical. When designing reliable real-world systems (e.g. cloud platforms and manufacturing lines), the faster the convergence, the faster the recovery after disturbances. When designing efficient sample-based algorithms (e.g. stochastic gradient descent (SGD) variants and MCMC), the faster the convergence, the faster the goal attainment. The rate of convergence also appears in MDP-related sample complexity results under the name "mixing time". Although convergence rate estimation is critically important, estimating the convergence rate of even a mildly complex chain can be extremely difficult.Over the last three decades, significant efforts have been made to bound the convergence of general state-space Markov chains. Most of these works utilize a pair of drift and minorization conditions (D&M) to bound the convergence in terms of the total variation (TV) distance (Meyn et al., 1994; Rosenthal, 1995; Jarner and Roberts, 2002; Douc et al., 2004; Baxendale, 2005; Andrieu et al., 2015) . The drift condition forces the chain to move towards a selected region. On such a region, the minorization condition allows the chain to regenerate or to couple with a stationary version of the chain. This analysis tends to produce overly conservative TV bounds, especially in high-dimensional settings; see (Qin and Hobert, 2021) for a discussion.The Wasserstein distance, as a measure of convergence to equilibrium, can exhibit better dimension dependence (Qin and Hobert, 2022b) . In addition, many Markov chains of interest (e.g. constant step-size SGD minimizing convex loss on finite datasets) converge in Wasserstein distance but not in TV distance. Consequently, bounding convergence in Wasserstein distance has steadily gained popularity over the years (Gibbs, 2004; Hairer et al., 2011; Butkovsky, 2014; Durmus and Moulines, 2015; Durmus et al., 2016; Qin and Hobert, 2022a) . Most of these works replace the minorization condition with a contraction condition (D&M becomes D&C). After returning to a selected region, two copies of the chain tend to become closer to each other. Both D&M and D&C enforce two conditions in two respective regions. However, partitioning the state space into two distinct regions often leads to suboptimal rates.Recently, (Qu et al., 2023) introduce the so-called contractive drift condition (CD), a single condition enforced on the entire state space, to explicitly bound the convergence in Wasserstein distance. A special case of CD dates back to (Steinsaltz, 1999) . By verifying CD, (Qu et al., 2023) establish parametrically sharp convergence bounds for stylized Markov chains arising from queueing theory and stochastic optimization (e.g. revealing how step-size, heavy-tailed gradient noise, growth rate and local curvature of objectives affect the convergence of stylized SGD). Although CD may generate better bounds than D&M and D&C for stylized chains (e.g. SGD with iid gradient noise), these methods are generally intended as theoretical tools that can provide closed-form convergence bounds for structured models. For more realistic, less structured chains, computational rather than analytical methods are needed. However, despite of the rapid development of computational power in the past decade, the convergence analysis of general state-space Markov chains is still in the pen-and-paper age.To launch computational Markov chain convergence analysis, we need a key to switch on the deep learning engine. This paper introduces the Deep Contractive Drift Calculator (DCDC) is the first general-purpose sample-based algorithm for bounding the convergence of general state-space Markov chains. There are two key ideas we develop. The first is to observe that CD, an inequality by definition, is actually an equality by nature (if the inequality has a solution, then the corresponding equality also has a solution). Thus, we introduce the Contractive Drift Equation (CDE), an integral equation the solution of which leads to an explicit convergence bound. For the second part, inspired by the success of physics-informed neural networks (PINNs) in solving PDEs (Sirignano and Spiliopoulos, 2018; Raissi et al., 2019) , we develop an efficient neural-network-based CDE solver. By combining these two components, DCDC solves CDEs by training neural networks and converts solutions into explicit convergence bounds. DCDC demonstrates the potential of computer-assisted convergence analysis and bridges the gap between deep learning and a traditionally challenging area of mathematical analysis.In high-dimensional spaces, PINNs minimize the integrated residual of a PDE via SGD to find a continuously differentiable function that approximately satisfies the PDE. When applying this idea to solve a CDE, an integral equation, the solving procedure becomes more natural in the following two ways. First, we only assume that the CDE solution is Lipschitz continuous, and neural networks are inherently Lipschitz continuous. Second, as SGD is already used to handle the integrated residual, we can simultaneously use it to handle the integral in the CDE. After approximately solving the CDE, DCDC needs to convert the solution into a convergence bound, which requires that the solution is uniformly accurate with high probability. This is different from PINNs in the PDE literature since the accuracy is mainly measured in the L 2 sense.The CDE solution is a new type of Lyapunov function that provides explicit convergence rates for random dynamical systems. For deterministic dynamical systems, traditional Lyapunov functions play central roles in establishing stability; see (Pukdeboon, 2011 ) for a review. There is a substantial literature on computing traditional Lyapunov functions via neural networks; see (Liu et al., 2023) and references therein. As pointed out in (Dawson et al., 2023) , a survey on certificate learning, learned (traditional) Lyapunov functions provide safety certificates for learned control policies (on deterministic dynamical systems). For the control of random dynamical systems, DCDC not only generates safety certificates (CDE solutions) but also quantifies safety levels (convergence rates). Control and performance evaluation of random dynamical systems have become a staple in contemporary data-driven decision making systems, thus underscoring the importance of DCDC.In short, we summarize our contributions as follows:• We introduce the Deep Contractive Drift Calculator (DCDC), the first general-purpose end-to-end approach that enables the use of deep learning to bound the convergence rate of general state-space Markov chains.• We perform sample complexity analysis and use DCDC to generate convergence bounds for realistic Markov chains arising in operations research as well as machine learning. • Our DCDC approach discovers features that are exploited by techniques developed to study CDs by closed-form methods, such as the wedge shape and the boundary removal technique discussed in (Qu et al., 2023) .Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?Answer: [Yes] Justification: All claims are justified.Guidelines:• The answer NA means that the abstract and introduction do not include the claims made in the paper. • The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers. • The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings. • It is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper. 
contractive drift equation: 
 Let X be a Markov chain on X ⊂ R d , with random mapping representationX n+1 = f n+1 (X n ), n = 0, 1, 2, . . .where f n 's are iid copies of f , a locally Lipschitz random mapping from X to itself (with probability one, f : X → X is locally Lipschitz).Example. Let α be a positive constant and Z be a square integrable random variable. The SGD with step-size α to solvemin x E(x -Z) 2 /2 is X n+1 = X n -α(X n -Z n+1 ) where Z n+1 's are iid copies of Z, so the corresponding random mapping is f (x) = x -α(x -Z).Understanding the long-term behavior of X requires estimating how fast X n converges to X ∞ (equilibrium) as n → ∞. The difference between the two distributions is quantified by either total variation (TV) distance or Wasserstein distance. Representative TV convergence bounds (e.g., T V (X n , X ∞ ) ≤ Cr n ) can be found in (Meyn et al., 1994; Rosenthal, 1995; Baxendale, 2005) .Representative Wasserstein convergence bounds (e.g., W (X n , X ∞ ) ≤ Cr n ) can be found in (Hairer et al., 2011; Durmus and Moulines, 2015; Qin and Hobert, 2022a) . These analytical methods can only handle stylized (structured) Markov chains. The goal of this paper is to introduce the first computational method that can handle realistic (less structured) Markov chains.The first step to achieve the goal is introducing the contractive drift equation (CDE). The local Lipschitz constant of f at x ∈ X is defined asDf (x) ∆ = lim δ→0 sup x ′ ,x ′′ ∈B δ (x) ∥f (x ′ ) -f (x ′′ )∥ ∥x ′ -x ′′ ∥where ∥•∥ is the Euclidean norm andB δ (x) = {x ′ : ∥x ′ -x∥ < δ}. If f is differentiable, then Df (x) = lim h→0 sup v:∥v∥=1 ∥f (x + hv) -f (x)∥ h = sup v:∥v∥=1 ∥∇f (x)v∥ = ∥∇f (x)∥where ∇f is the Jacobian matrix of f and ∥•∥ becomes the spectral norm when applying to matrices. Basically, Df (x) describes how expansive or contractive f is around x. With these notations, the contractive drift condition (CD) in (Qu et al., 2023) that leads to computable convergence bounds isKV (x) ∆ = EDf (x)V (f (x)) ≤ V (x) -U (x), x ∈ X (1)where V, U : X → R + are bounded away from zero. In the rest of this paper, we adopt the convention that all functions denoted by U are positive and bounded away from zero, i.e. inf U > 0. We use E x to denote the expectation operator conditional on X 0 = x. In (1), the subscript is omitted as the initial location is clear. By replacing "≤" with "=" in (1), the contractive drift equation (CDE) is KV = V -U , for which we establish the following existence and uniqueness results. All proofs are in the appendix. Theorem 1. Fix U and suppose that KW ≤ W -U has a non-negative finite solution W * . ThenV * (x) ∆ = E x ∞ k=0 U (X k ) k l=1 Df l (X l-1 ) , x ∈ X (2)is finite and satisfies KV * = V * -U . Furthermore, KV = V -U has at most one bounded solution.Remark. This V * can be interpreted as an average space-discounted cumulative reward. Imagine a swarm of agents moving according to f . For an agent at x, if Df (x) < 1 (contraction), then after f is applied, there will be more agents around this agent. If all agents around f (x) share a total reward U (f (x)), then the reward for each of them is discounted. From the perspective of a particular agent, the procedure is like collecting reward within a shrinking ball.3 Deep Contractive Drift CalculatorQuestion: Does the paper discuss the limitations of the work performed by the authors?Answer: [Yes] Justification: We have a "Limitations" section.Guidelines:• The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper. • The authors are encouraged to create a separate "Limitations" section in their paper.• The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be. • The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with a few runs. In general, empirical results often depend on implicit assumptions, which should be articulated. • The authors should reflect on the factors that influence the performance of the approach.For example, a facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or a speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon. • The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size. • If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness. • While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an important role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations. 
why do we introduce cde?: 
 Physics-informed neural networks (PINNs) solve a PDE by minimizing its integrated residual (Sirignano and Spiliopoulos, 2018; Raissi et al., 2019) . If we want to use this idea to solve KV ≤ V -U , then the integrated residual isl(θ) ∆ = X (KV θ (x) -V θ (x) + U (x)) + h(x)dxwhere h is a positive density and {V θ : θ ∈ Θ} is a neural network. Note that the residual at x is positive if and only ifKV θ (x) > V θ (x) -U (x). By letting X 0 have distribution h, l(θ) = E [E [Df 1 (X 0 )V θ (f 1 (X 0 )) -V θ (X 0 ) + U (X 0 )|X 0 ]] + ,which is an expectation of a non-linear function of a conditional expectation. Minimizing l(θ) is a conditional stochastic optimization problem (CSO). In CSO, the sample-average gradient is biased (Hu et al., 2020b) , which leads to a high sample complexity for convergence (Hu et al., 2020a) .Fortunately, if we aim at solving KV = V -U (CDE) instead of KV ≤ V -U (CD), then there exists a simple unbiased gradient estimator. Now we briefly derive this estimator. For a CDE, the integrated residual becomesl(θ) ∆ = X (KV θ (x) -V θ (x) + U (x)) 2 h(x)dx =E [E [Df 1 (X 0 )V θ (f 1 (X 0 )) -V θ (X 0 ) + U (X 0 )|X 0 ]] 2 with its gradient l ′ (θ) =2E [E [Df 1 (X 0 )V θ (f 1 (X 0 )) -V θ (X 0 ) + U (X 0 )|X 0 ] E [Df 1 (X 0 )V ′ θ (f 1 (X 0 )) -V ′ θ (X 0 )|X 0 ]] =2EE [[Df 1 (X 0 )V θ (f 1 (X 0 )) -V θ (X 0 ) + U (X 0 )] [Df -1 (X 0 )V ′ θ (f -1 (X 0 )) -V ′ θ (X 0 )] |X 0 ] =2E [[Df 1 (X 0 )V θ (f 1 (X 0 )) -V θ (X 0 ) + U (X 0 )] [Df -1 (X 0 )V ′ θ (f -1 (X 0 )) -V ′ θ (X 0 )]]where f 1 and f -1 are iid copies of f while V ′ θ = dV θ /dθ is computed via backpropagation. This expression allows us to estimate l ′ (θ) without any bias. In summary, the inequality (CD) is enough to bound the convergence, but the equality (CDE) turns out to be easier to establish (via deep learning).Given the above discussion, a standard application of SGD is enough to simultaneously handle the integrated residual as well as the integral in the CDE, resulting in the following simple algorithm, Deep Contractive Drift Calculator, the first general-purpose sample-based algorithm to bound the convergence of general state-space Markov chains.Algorithm 1 Deep Contractive Drift Calculator (DCDC) Require: Step-size α, number of iterations T , neural network {V θ : θ ∈ Θ}, initialization θ 0 for t ∈ {0, ..., T -1} do sample (X 0 , f 1 , f -1 ) compute l′ (θ t ) as 2 [Df 1 (X 0 )V θt (f 1 (X 0 )) -V θt (X 0 ) + U (X 0 )] Df -1 (X 0 )V ′ θt (f -1 (X 0 )) -V ′ θt (X 0 )update θ t+1 = θ t -α l′ (θ t ) (SGD or its variants) end for convert V θ T into a convergence bound (Theorem 3 and Theorem 4)The conversion will be discussed in the next two subsections. In the current subsection, we show the validity of approximating CDE solutions via neural networks. In the following, we use ∥•∥ ∞ to denote the sup norm of functions on X .Theorem 2. If X is compact, ∥EDf ∥ ∞ is finite, and V * in (2) is finite and continuous, then for any ϵ > 0, there exists a neural network {V θ : θ ∈ Θ} and its realization V θ * such that∥KV θ * -V θ * + U ∥ ∞ < ϵ.Although DCDC solves CDEs on compact sets, it can be applied to Markov chains on non-compact sets that have compact absorbing sets (e.g. SGD for regularized problems). For chains without a compact absorbing set, extending DCDC to bound their convergence is left for future research, but here we describe a natural strategy to do so. In general, a Markov chain spends most of its time on some large compact set C where the chain may have complex dynamics. When the chain is outside C, it typically has a strong tendency to return. Therefore, to extend DCDC, we can (i) search some parametric family (e.g. V A (x) = x ⊤ Ax) to establish a CD outside C (capturing the return tendency);(ii) apply DCDC to obtain a CDE solution on C (capturing the complex dynamics); (iii) stitch them together to obtain a global CD. Comparing the large set here with the small set (Meyn and Tweedie, 2009) in D&M or D&C illustrates the advantage of computational methods over analytical ones. The size of the large set is determined by the approximation capability of neural networks, but the size of the small set is determined by the minorization or contraction condition (the two conditions often require the small set to be very small).Now we discuss how to convert KV ≤ V -U into convergence bounds with exponential rates in Wasserstein distance. To begin, we recall the definition of the Wasserstein distance. Let P(X ) be the set of probability measures on X equipped with its Borel sigma-algebra. The Wasserstein distance between µ, ν ∈ P(X ) isW (µ, ν) ∆ = inf π∈C(µ,ν) X ×X ∥x -y∥ π(dx, dy) where C(µ, ν) ∆ = {π ∈ P(X × X ) : π(•, X ) = µ(•), π(X , •) = ν(•)} is the set of all couplings of µ and ν. Given two random variables Z 1 and Z 2 , we use W (Z 1 , Z 2 ) to denote the Wasserstein distance between their marginal distributions. Theorem 3. Suppose that X is convex and thatKV ≤ V -U holds with sup V < ∞. If E ∥X 0 -X 1 ∥ < ∞, then X has a unique stationary distribution X ∞ with W (X n , X ∞ ) ≤ Cr n , r ∆ = 1 -inf U/ sup V, C ∆ = E ∥X 0 -X 1 ∥ V (X 0 + Ũ (X 1 -X 0 )) inf U • (inf V / sup V )where Ũ is a U [0, 1] random variable independent of X 0 and X 1 .Given U , the exponential rate r is determined by the magnitude of V . The smaller the V , the faster the convergence. Given X 0 , the pre-multiplier C can be easily computed by simulating the first transition (from X 0 to X 1 ).In Theorem 3 of (Qu et al., 2023) , convergence bounds with exponential rates are straightforwardly derived from KV ≤ rV where r < 1, so one might wonder why we need the less straightforward Theorem 3 here. This is because KV ≤ rV is not suitable for PINN-like solvers. In Theorem 3, we solve KV = V -U and compute the exponential rate r from the solution V . However, for KV = rV , we need the answer (the exponential rate r) to write down the question (the equation to solve and the corresponding loss to minimize), which is circular. Of course, we may try solving KV = rV for different values of r, but it turns out that it is very hard for DCDC to converge even for very conservative (close to 1) r's. Here is an explanation. Unlike KV = V -U , which has a solution as long as KV ≤ V -U has one (Theorem 1), KV = rV may not have a solution even when KV ≤ rV has one. However, it is not hard to show that KV = rV -r has a (formal) solutionV r (x) ∆ = E x ∞ k=0 (1/r) k k l=1 Df l (X l-1 ) , x ∈ X .Comparing with V * in (2), U (X k ) is replaced by exponentially exploding (1/r) k . Back to KV = rV , its solution (if there is any) should be the above expression without the summation but with k → ∞ (as a limit), which suggests that the solution may have a large magnitude, making it difficult to approximate.Now we discuss how to generate convergence bounds with polynomial rates using DCDC. The key is to iteratively solve a sequence of CDEs. For example, given V 0 , we first solveKV 1 = V 1 -V 0 to obtain V 1 . Then we solve KV 2 = V 2 -V 1 to obtain V 2 .These two CDEs together lead to an O(1/n) convergence bound. Theorem 4. Suppose that X is convex and that there exist positive functionsV 0 , V 1 , . . . , V m such that 0 < inf V 0 < sup V m < ∞ and KV k+1 ≤ V k+1 -V k for k = 0, . . . , m -1. If E ∥X 0 -X 1 ∥ < ∞, then X has a unique stationary distribution X ∞ with W (X n , X ∞ ) ≤ E ∥X 0 -X 1 ∥ V m (X 0 + Ũ (X 1 -X 0 )) inf V 0 • m-1 k=1 (1 + n/k) where Ũ is a U [0, 1] random variable independent of X 0 and X 1 .The expectation in the numerator can be easily computed by simulating the first transition, while the product in the denominator is basically n m-1 as n → ∞.In Theorem 1 of (Qu et al., 2023) , convergence bounds with polynomial rates (O(1/n m-1 )) are derived from KV ≤ V -U 1/m V 1-1/m paired with KU ≤ U , so one might wonder why we need so many CDs in Theorem 4 here. This is because KV ≤ V -U 1/m V 1-1/m is designed for the pen-and-paper setting where directly establishing a sequence of CDs is difficult. GivenKV ≤ V -U 1/m V 1-1/m ,many inequalities are applied to extract a CD sequence from this single special CD, resulting in large constants in convergence bounds. DCDC makes it possible to directly establish a sequence of CDs (by consecutively solving CDEs). In this setting, we can use Theorem 4 to obtain better convergence bounds. To be specific, compared with our Theorem 4, the result in (Qu et al., 2023) has an extra factor m m /m!.Question: For each theoretical result, does the paper provide the full set of assumptions and a complete (and correct) proof?Answer: [Yes] Answer: [Yes] Justification: Code with instructions is included in the submission.Guidelines:• The answer NA means that paper does not include experiments requiring code. • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
sample complexity: 
 As a numerical solver, DCDC solves CDEs approximately. Let Ṽ = V θ T be the output of DCDC. We should not expect K Ṽ = Ṽ -U to hold exactly. Even if Ṽ is an exact solution, the exactness is hard to verify as K Ṽ is an expectation and the domain X ⊂ R d is not a finite set. As establishing convergence bounds requires CDs to exactly hold everywhere, given N iid copies of f to estimate K and M = {x 1 , . . . , x M } uniformly sampled from X , we can (i) establishKN Ṽ (x) ∆ = 1 N N k=1 Df k (x) Ṽ (f k (x)) ≤ Ṽ (x) -Ũ (x), x ∈ Mwhere Ũ may be smaller thanU (e.g. if Ṽ is supposed to solve Ṽ -K Ṽ = U ≡ 1, then Ũ ≡ inf M [ Ṽ -KN Ṽ ]); (ii) claim that K Ṽ ≤ Ṽ -Ũ + ϵ holds everywhere with probability at least 1 -δ where ϵ, δ > 0; (iii) convert K Ṽ ≤ Ṽ -Ũ + ϵ into a convergence bound.To have M, N large enough to make the claim in (ii), we need the following sample complexity result.Theorem 5. Suppose that (i) X is compact; (ii) V, U are bounded and Lipschitz; (iii) EDf 2 + ED 2 f < ∞ where Df is the Lipschitz constant of f and D 2 f is the Lipschitz constant of Df . Given ϵ, δ > 0, we can choose M = O(log(1/ϵ)/(δϵ d )) and N = O(1/(δϵ 2 )) to have P sup x∈X [KV (x) -V (x) + U (x)] ≤ sup x∈M KN V (x) -V (x) + U (x) + ϵ > 1 -δ.Since the exponential rate of convergence in Theorem 3 is r = 1 -inf U/ inf V , Theorem 5 also provides the sample complexity for estimating the exponential rate. Specifically, with probability at least 1 -δ, the exponential rate rM,N computed from KN V ≤ V -U on M, which may not be a valid exponential rate, is ϵ-close to a valid exponential rate r * (given by KV ≤ V -U + ϵ on X ).It is worth noting that in terms of sample complexity, Theorem 5 guarantees a DCDC certificate (and thus a convergence bound to stationarity) with high probability with an efficient parametric O(1/N 1/2 ) rate in terms of the number of samples (namely, the bound holds with high probability up to an error of order O(1/N 1/2 )). Once samples are generated, M = Õ(1/ϵ d ) points are chosen for the empirical evaluation. Thus, the total complexity (both in terms of number of evaluations and number of samples is O(1/ϵ 2 ) + Õ(1/ϵ d ). A related literature on parametric integration (i.e. learning a Markov transition kernel that maps Lipschitz functions to continuous functions on the d-dimensional cube) provides a lower bound of order Õ(1/ϵ d ), (Heinrich and Sindambiwe, 1999) . Although these results are suggestive, they cannot be applied directly because we assume a random mapping representation, which provides additional structure. We plan to study the lower bounds in future work.5 Numerical Examples 
mini-batch sgd for logistic regression with regularization: 
 Having established the theoretical foundation of DCDC, we now utilize it to generate convergence bounds for Markov chains of interest that are too hard for pen-and-paper analysis. To begin, we bound the convergence of a constant step-size mini-batch SGD that minimizes the cross-entropy loss over a finite dataset with L 2 regularization.Let (x 1 , y 1 ),. . . ,(x m , y m ) be m data points where x i ∈ [-1/2, 1/2] 2 and y i ∈ {0, 1}. To perform regularized logistic regression, we want to choose b ∈ R 2 to minimize- 1 m m i=1 (y i log p i + (1 -y i ) log(1 -p i )) + λ 2m ∥b∥ 2 , p i = σ(b ⊤ x i ) = 1 1 + exp(-b ⊤ x i )where λ > 0 is the regularization parameter. The random mapping representation of the corresponding SGD with step-size α and batch-size β isf (b) =b(1 -λα/m) + (α/β) i∈B y i -σ(b ⊤ x i ) x iwhere B with |B| = β is uniformly sampled from {1, . . . , m}. Thanks to the regularization, the chain has a compact absorbing set. In fact, the chain cannot escape from B m/(λ√ 2) (0). The local Lipschitz constant of f (b) is Df (b) = (1 -λα/m)I -(α/β) i∈B σ ′ (b ⊤ x i )x i x ⊤ i ≤ 1 -λα/mwhere ∥A∥ = sup v:∥v∥=1 ∥Av∥ is the spectral norm. This demonstrates that the L 2 regularization makes the chain contractive ∥f (b 1 ) -f (b 2 )∥ ≤ (1 -λα/m) ∥b 1 -b 2 ∥. However, since the regularization parameter is chosen via cross-validation in a separate validation process, it is useful to obtain a contraction rate that is uniform in the regularization parameter. This rate is brought by the second term in Df (b) -we refer to this contribution as the intrinsic convergence rate. However, it is challenging to analyze the spectrum of this state-dependent data-based random matrix, so we need DCDC. The code is available in the supplementary material. Each training procedure in this paper was completed within ten minutes on an M2 MacBook Air with 8GB RAM.For the dataset, we set m = 100 and uniformly generate 100 x i 's. For each x i , its label y i follows Ber(0.9) or Ber(0.1), depending upon which coordinate of x i is larger. For the SGD, we set the regularization parameter λ = 1, step-size α = 0.1, and batch-size β = 10. For DCDC, we run 1M Adam steps to train a single-layer network with width 1000 and sigmoid activation. We also experiment with deeper networks with the same amount of neurons, and the results are similar.As demonstrated in Figure 5 .1, the single-layer network can already accurately solve the CDE KV = V -0.1. The learned solution Ṽ is on the left while the estimated difference K Ṽ -Ṽ is on the right. Aiming at KV ≤ V -0.1, we get K Ṽ ≤ Ṽ -0.0986. This leads to exponential rate 1 -1.07 × 10 -3 (Theorem 3) where 1 × 10 -3 corresponds to the regularization contribution, while 7 × 10 -5 corresponds to the intrinsic rate. Now we briefly discuss how the surface in Figure 5 .1 (left) leads to the intrinsic rate. From the expression of Df (b), we know that the intrinsic contraction concentrates around the center. To make it contribute to the overall convergence, it needs to be spread. The surface in Figure 5 .1 (left) provides the media to spread: (i) for points not at the center, the sunken surface creates a drift E x V (X 1 ) < V (x); (ii) however, for points at the center, the sunken surface creates an anti-drift E x V (X 1 ) > V (x), but it is overcome by the strong contraction E x Df 1 (x)V (X 1 ) < V (x). In this way, the strong contraction is spread (in the form of drift) to overall improve the contractive drift, which leads to the intrinsic rate. To conclude this example, when X 0 = 0, we compute the pre-multiplier C = 8.1, which leads to convergence bound W (X n , X ∞ ) ≤ 8.1(1 -1.07 × 10 -3 ) n .In the above SGD example, contraction plays the leading role. Now we consider a tandem fluid network (Kella and Whitt, 1992) where drift plays the leading role. Let s 1 and s 2 be two stations with buffer capacity c that can process fluid workload at rates r 1 and r 2 , respectively. External fluid only arrives at s 1 and is processed by s 1 then s 2 . Assume that the external input follows a compound renewal process where a random amount of fluid Z arrives after a random length of time T has passed since the last arrival. If r 1 < r 2 , then s 2 is always empty, so we let r 1 > r 2 . Let X be the remaining workload vector after each arrival. Its random mapping representation isf (x 1 , x 2 ) = (min((x 1 -r 1 T ) + +Z, c), (min(x 2 +(r 1 -r 2 ) min(T, x 1 /r 1 ), c)-r 2 (T -x 1 /r 1 ) + ) + )where x 1 decreases at rate r 1 until it is empty while x 2 increases at rate (r 1 -r 2 ) until x 1 is empty. Basically, within [0, c] 2 , the chain follows a northwest-then-south path for time T and then jumps east by amount Z. This chain has simple local Lipschitz constant Df (x 1 , x 2 ) = I(T ≤ (x 1 + x 2 )/r 2 ), obtained as an infinitesimal ball around (x 1 , x 2 ) collapses to a single point when the system is depleted before the next arrival. As a result, drift plays the leading role as contraction only happens around the origin.For the tandem network, we set the buffer capacity c = 1, processing rates (r 1 , r 2 ) = (1.1, 1.0), interarrival time T ∼ U [0, 0.2] and arriving amount Z ∼ U [0, 0.1] (the stability condition is EZ < r 2 ET ). For DCDC, we run 1M Adam steps to train a double-layer network with width 40 and sigmoid activation.Although a slightly deeper network is trained, the result in Figure 5 .2 (left) is almost a plane. Now we briefly explain why this is the correct solution. First, note that as long as the stability condition EZ < r 2 ET holds, the total workload V (x 1 , x 2 ) = x 1 + x 2 is the most natural Lyapunov function such that E x V (X 1 ) -V (x) = E(-r 2 T + Z) < 0 holds when x is far away from the boundary. Second, the "boundary removal technique" introduced in (Qu et al., 2023) shows that the above drift can be extended to the boundary as a contractive drift E x Df (x) V (X 1 ) -V (x) = E(-r 2 T + Z) < 0 as if the boundary (that causes anti-drift) never exists. The plane in Figure 5 .2 (left) demonstrates that DCDC has already mastered the above two steps! Again, we conclude this example with convergence bound W (X n , X ∞ ) ≤ 5.67(1 -0.017) n when X 0 = 0. Lyapunov functions are usually denoted by V in the literature, and V typically represents the shape of Lyapunov functions. As mentioned in the introduction, the CDE solution is a new type of Lyapunov function. In most cases, it is also V -shaped, representing the drift towards some contractive region. However, Markov chains sometimes exhibit neither drift nor contraction, such as when SGD is stuck in a non-strongly-convex basin or when the water level of the Moran dam (Stadje, 1993) is neither too low nor too high. Here, we use the simplest example, a two-sided regulated random walk f (x) = max(min(x + Z, 1/2), -1/2) with Z ∼ U [-1/3, 1/3], to illustrate that DCDC discovers upside-down V -shaped Lyapunov functions to address the above issue. In [-1/6, 1/6], the chain exhibits neither drift (Z is symmetric) nor contraction (±1/2 boundaries are not reachable in one step). The wedge in Figure 5 .3 (left) creates an artificial drift to maintain the CD.In (Qu et al., 2023) , a similar function is introduced as a tool to study stylized non-strongly-convex SGD. DCDC not only discovers this tool but also makes the wedge meaningful. As mentioned in the remark below Theorem 1, the CDE solution generated by DCDC represents an average spacediscounted cumulative reward, where an agent collects reward within a shrinking ball. Why does starting from the middle lead to the highest reward? Because the ball starting there has the longest lifespan before hitting the boundary and collapsing into a single point.Finally, to demonstrate the potential of DCDC to accurately recover exact convergence rates, we examine a class of d-dimensional autoregressive processesf (x) = Hx + Z, x, H, Z ≥ 0 (3)where random vector Z is integrable and constant matrix H is symmetric with all its eigenvalues in (0, 1). The exact convergence rate of this Markov chain in Wasserstein distance can be explicitly computed by pen and paper. Proposition 1. Let X be the Markov chain defined by (3). If X 0 = 0, thenE ∥H n Y ∥ / √ d ≤ W (X n , X ∞ ) ≤ E ∥H n Y ∥ , Y = ∞ k=1 H k-1 Z kwhere Z k 's are iid copies of Z. Let λ be the largest eigenvalue of H. Then E ∥H n Y ∥ = Θ(λ n ) as long as Y does not concentrate on the orthogonal complement of the eigenspace associated with λ.For the autoregressive process, let d = 3 and H = 0.4 0.2 0.1 0.2 0.5 0.2 0.1 0.2 0.6 with λ = 0.850. Let Z be uniformly sampled from B 1 (0) ∩ R 3 + . Note that the resulting Markov chain cannot escape from B 10 (0) ∩ R 3 + . After plugging into the simulator of (3), DCDC generates V ≡ 0.668, which recovers the exact convergence rateKV = V -U = (1 -U/V )V = (1 -0.1/0.668)V = 0.850V = λV. 
conclusions: 
 We introduce DCDC, a potent framework that enables the use of deep learning techniques to tackle the problem of estimating convergence to stationarity of complex, general state-space Markov chains. Our approach unlocks the key to using scalable data-driven tools to tackle this important problem. In future work, we plan to use these results in the context of general state-space reinforcement learning, control of ergodic systems, and related applications by employing the CD condition as a policy regularizer.Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?Answer: [Yes] Justification: We present enough details in the paper and full details in the code.Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
limitations: 
 DCDC solves CDEs on compact spaces. A potential strategy to handle non-compact spaces is discussed in Section 3.2. Sample complexity lower/upper bounds are not studied in this paper and they are left for future research.Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?Answer: [Yes] Justification: We report the average error and its standard deviation in the numerical experiments.Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes] Justification: As mentioned in the paper, an M2 MacBook Air suffices with 8GB RAM suffices.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines?Answer: [Yes] Justification: Our research conforms with the Code of Ethics.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed? 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: The paper poses no such risks Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [NA] Justification: The paper does not use existing assets Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. 
 
 Appendix: Proof of Theorem 1. Note thatW * (x) ≥KW * (x) + U (x) =EDf 1 (x)W * (f 1 (x)) + U (x) ≥EDf 1 (x)(KW * (f 1 (x)) + U (f 1 (x))) + U (x) =EDf 1 (x)E [Df 2 (f 1 (x))W * (f 2 (f 1 (x)))|f 1 ] + EDf 1 (x)U (f 1 (x)) + U (x) =EDf 1 (x)Df 2 (f 1 (x))W * (f 2 (f 1 (x))) + EDf 1 (x)U (f 1 (x)) + U (x) • • • ≥E x W * (X n ) n k=1 Df k (X k-1 ) + n-1 k=0 E x U (X k ) k l=1 Df l (X l-1 ) ≥ n-1 k=0 E x U (X k ) k l=1Df l (X l-1 ) .As n → ∞, we have V * ≤ W * < ∞ andKV * (x) = EDf 1 (x)E X1 ∞ k=0 U (X k+1 ) k l=1 Df l+1 (X l ) = E ∞ k=0 U (X k+1 ) k l=0 Df l+1 (X l ) = E ∞ k=1 U (X k ) k l=1 Df l (X l-1 ) = V * (x) -U (x).Let V * be another solution of KV = V -U . Similar to W * ,V * (x) = E x V * (X n ) n k=1 Df k (X k-1 ) + n-1 k=0 E x U (X k ) k l=1Df l (X l-1 ) .(4)As n → ∞, we have V * ≥ V * . If V * , and hence V * , is unbounded, then KV = V -U doesn't have bounded solution. If V * , and hence V * , is bounded, we claim that they are the same solution. It suffices to show that the first term on the RHS of (4) vanishes as n → ∞. This is true becauseV * (x) < ∞ ⇒ E x U (X n ) n k=1 Df k (X k-1 ) → 0 ⇒ E x V * (X n ) n k=1 Df k (X k-1 ) → 0where the last step follows from sup V * < ∞ and inf U > 0.Proof of Theorem 2. Since continuous functions on compacts sets are bounded, by Theorem 1, V * is the unique continuous solution of KV = V -U. By the universal approximation theorem (Cybenko, 1989) , there exists a single-layer neural network with sigmoid activation {V θ : θ ∈ Θ} (Θ is some Euclidean space) and its realizationV θ * such that ∥V θ * -V * ∥ ∞ < ϵ/( L + 1) where L = ∥EDf ∥ ∞ . Then ∥KV θ * -V θ * + U ∥ ∞ ≤ ∥KV * -V * + U ∥ ∞ + ∥V θ * -V * ∥ ∞ + ∥KV θ * -KV * ∥ ∞ <0 + ϵ/( L + 1) + Lϵ/( L + 1) =ϵ.Proof of Theorem 3. The setting introduced in Section 2 is a special case of the setting in (Qu et al., 2023) , allowing us to directly invoke the results from that work. In particular, as X is assumed to be convex, the intrinsic metric used in (Qu et al., 2023) reduces to the Euclidean metric. From KV ≤ V -U , we haveKV ≤ V -U ≤ V -U • V / sup V ≤ rV, r = 1 -inf U/ sup V.By Theorem 3 in (Qu et al., 2023) ,W (X n , X ∞ ) ≤ 1 inf V r n 1 -r E ∥X 0 -X 1 ∥ V (X 0 + Ũ (X 1 -X 0 )) = Cr n where C = E ∥X 0 -X 1 ∥ V (X 0 + Ũ (X 1 -X 0 )) inf U • (inf V / sup V )and Ũ is a U [0, 1] random variable independent of X 0 and X 1 .Proof of Theorem 4. The proof is similar to the proof of Theorem 1 in (Qu et al., 2023) , but in our specific setting, the proof becomes much simpler notation-wise. As in the proof of our Theorem 1, we haveV 1 (x) ≥ ∞ n=0 K n V 0 (x), K n V 0 (x) = E x V 0 (X n ) n l=1Df l (X l-1 ) .Following the same induction process as in (Qu et al., 2023) , we haveV m (x) ≥ ∞ n=0 c m,n K n V 0 (x), c m,n = n + m -1 m -1 . Let F n = f n • • • • • f 1 and Fn = f 1 • • • • • f nbe the n-fold forward and backward composition, respectively. Given f independent of anything else,f (x) x V m (y)dy ≥ ∞ n=0 c m,n f (x) x K n V 0 (y)dy = ∞ n=0 c m,n E f (x) x V 0 (F n (y)) n l=1 Df l (F l-1 (y)) dy ≥ ∞ n=0 c m,n E Fn(f (x)) Fn(x) V 0 (z)dz ≥ ∞ n=0 c m,n E Fn+1(x) Fn(x) V 0 (z)dzFor a particular n,f (x) x V m (y)dy ≥ ∞ n=n c m,n E Fn+1(x) Fn(x) V 0 (z)dz ≥c m,n E F∞(x) Fn(x) V 0 (z)dz ≥c m,n inf V 0 • E Fn (x) -F∞ (x) .By integrating with respect to the initial distribution X 0 ,E ∥X 0 -X 1 ∥ V m (X 0 + Ũ (X 1 -X 0 )) ≥c m,n inf V 0 • E Fn (X 0 ) -F∞ (X 0 ) ≥c m,n inf V 0 • W (X n, X ∞ ) = (n + m -1) . . . (n + 1) (m -1) . . . 1 inf V 0 • W (X n, X ∞ ) = m-1 k=1 n k + 1 inf V 0 • W (X n, X ∞ ).Proof of Theorem 5. As X is compact, without loss of generality, let X = [0, 1] d . The main goal is to find M, N such thatP sup x∈X [KV (x) -V (x) + U (x)] > sup x∈M KN V (x) -V (x) + U (x) + ϵ ≤ δ.This probability is bounded by the sum ofEQUATIONandEQUATIONTo bound (5), we need to bound the Lipschitz constant of KV -V + U . For x, y ∈ X ,|KV (x) -KV (y)| = |EDf (x)V (f (x)) -EDf (y)V (f (y))| ≤ |EDf (x)V (f (x)) -EDf (x)V (f (y))| + |EDf (x)V (f (y)) -EDf (y)V (f (y))| ≤EDf (x) |V (f (x)) -V (f (y))| + E |Df (x) -Df (y)| V (f (y)) ≤DV • EDf 2 • ∥x -y∥ + sup V • ED 2 f • ∥x -y∥ where DV is the Lipschitz constant of V . Then, the Lipschitz constant of KV -V + U is bounded by L ∆ = DV • EDf 2 + sup V • ED 2 f + DV + DU.Then ( 5) is bounded byP [0, 1] d ̸ ⊂ ∪ x∈M B r where r = ϵ/(2 L).To bound this probability, we divide the unit cube into C = (√ d/r) d = (2 L√ d/ϵ) d sub-cubes with edge length r/ √ d. Then [0, 1] d ̸ ⊂ ∪ x∈M B rimplies that there exists at least one sub-cube that does not contain any element of M. This is equivalent to failing to collect C different coupons within M draws. It is well-known that we need Θ( C log C) on average to collect C different coupons. By Markov inequality, we can choose M = O( C log C/δ) = O(log(1/ϵ)/(δϵ d )) to reduce the failure probability below δ/2.To bound (6), let x M = arg max x∈M [KV (x) -V (x) + U (x)]. By Chebyshev inequality,P KV (x M ) -V (x M ) + U (x M ) > sup x∈M KN V (x) -V (x) + U (x) + ϵ/2 M ≤P KV (x M ) -V (x M ) + U (x M ) > KN V (x M ) -V (x M ) + U (x M ) + ϵ/2 M =P KV (x M ) > KN V (x M ) + ϵ/2 M ≤ Var [Df (x M )V (f (x M ))|M] /N ϵ 2 /4 ≤ 4 sup V 2 • EDf 2 N ϵ 2 ≤δ/2 when N ≥ (8 sup V 2 • EDf 2 )/(δϵ 2 ) = O(1/(δϵ 2 )).Proof of Proposition 1. By X 0 = 0 and (3), we haveX n = n k=1 H n-k Z k and Xn = n k=1 H k-1 Z kwhere Xn is the backward chain. By definition,W (X n , X ∞ ) ≤ E X∞ -Xn = E ∞ k=n+1 H k-1 Z k = E H n ∞ k=1 H k-1 Z k = E ∥H n Y ∥ .Let Lip(1) be the family of functions on R d that are 1-Lipschitz with respect to the 2-norm ∥•∥. Let ∥•∥ 1 be the 1-norm. By the Kantorovich-Rubinstein duality,W (X n , X ∞ ) = sup h∈Lip(1) |Eh(X ∞ ) -Eh(X n )| ≥E X∞ 1 -Xn 1 / √ d =E X∞ -Xn 1 / √ d ≥E X∞ -Xn / √ d =E ∥H n Y ∥ / √ d,where the second line is because∥•∥ 1 / √ d is 1-Lipschitz with respect to ∥•∥, the third line is because X∞ ≥ Xn ≥ 0, and the fourth line is because ∥•∥ 1 ≥ ∥•∥.Let Λ be the diagonal matrix containing all the eigenvalues of H. Let Q be the orthogonal matrix containing all the (column) eigenvectors of H. ThenH n = QΛ n Q ⊤ = d i=1 λ n i q i q ⊤ i .Recall that λ is the largest eigenvalue of H. ThenE ∥H n Y ∥ /λ n = E d i=1 (λ i /λ) n q i q ⊤ i Y → E i∈I q i q ⊤ i Y , as n → ∞where {q i : i ∈ I} are the eigenvectors corresponding to λ. Since Z is integrable and Y does not concentrate on the orthogonal complement of the eigenspace associated with λ, the above limit is finite and positive.NeurIPS Paper ChecklistJustification: There is no societal impact of the work performed.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.• The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 


Decision: This paper proposes a general-purpose deep learning-based algorithm to compute the convergence rate of Markov chains to their stationary distribution. The proposed method consists of two components:

1. **Theoretical Tool:**  
   The first component is a theoretical tool that demonstrates how the convergence rate of the Markov chain can be exactly determined by solving a Contractive Drift Equation (CDE), derived from the contractive drift condition of the Markov chain.

2. **Network-Based Solver:**  
   The second component is a network-based solver for the CDE equation. The authors provide examples on Stochastic Gradient Descent (SGD) and other practical Markov chains to showcase how the convergence rate can be obtained for complex, real-world chains.

All reviewers agree that this work is well-written and presents a novel solution that bridges the gap between the theoretical analysis of Markov chains and data-driven algorithms. As suggested by the reviewers, it would be beneficial to include refined experiments that verify the convergence rate by running sufficient steps of the Markov chains.

Overall, this paper is a clear acceptance, and the idea of using deep learning to estimate the convergence of Markov chains deserves more attention from the community.

Reviews: review: paper_summary: paper_summary: #### Summary
The paper introduces the Deep Contractive Drift Calculator (DCDC), a novel sample-based algorithm for bounding the convergence rates of general state-space Markov chains to stationarity in Wasserstein distance. The method leverages deep learning to solve the Contractive Drift Equation (CDE), providing explicit convergence bounds. The paper includes theoretical analysis, sample complexity, and empirical validation on realistic Markov chains.
strength_weakness: strength_weakness: #### Weaknesses
1. **Computational Complexity**: The approach may be computationally intensive, particularly for high-dimensional state spaces. More discussion on computational efficiency and scalability would be beneficial.
2. **Comparison with Existing Methods**: While the paper discusses theoretical advantages, empirical comparisons with existing state-of-the-art methods for convergence analysis are limited.
3. **Generality**: The method is demonstrated on specific types of Markov chains. Extending the empirical validation to a broader range of applications would strengthen the paper.
4. **Sample Complexity**: Although the sample complexity is analyzed, practical guidelines for choosing sample sizes in different scenarios would enhance the utility of the method.
questions: questions: #### Questions
1. **Computational Complexity**:
    - Could you provide more details on the computational complexity of the DCDC algorithm, particularly for high-dimensional state spaces? How does the method scale with increasing dimensions?

2. **Comparison with Existing Methods**:
    - How does the DCDC method compare empirically with existing state-of-the-art methods for bounding convergence rates of Markov chains? Are there specific scenarios where DCDC significantly outperforms these methods?

3. **Generality**:
    - The method is validated on specific types of Markov chains. Do you foresee any challenges in applying DCDC to other types of Markov chains, such as those with more complex dynamics or in higher-dimensional spaces?

4. **Sample Complexity**:
    - While you provide a theoretical analysis of sample complexity, can you offer practical guidelines or heuristics for choosing sample sizes in different applications? How sensitive is the method to the choice of sample size?

5. **Practical Applications**:
    - Can you discuss potential practical applications of the DCDC method in more detail? For instance, how might this method be applied in real-world scenarios such as reinforcement learning or stochastic optimization?

6. **Assumptions and Limitations**:
    - The paper discusses some assumptions and limitations. Could you elaborate on the key assumptions that are critical for the theoretical results, and how robust the method is to violations of these assumptions?
limitations: limitations: no
paper_summary: paper_summary: The paper studies the problem of convergence rate analysis for general state-space Markov chains. They propose Deep Contractive Drift Calculator (DCDC), the first general-purpose sample-based algorithm for bounding the convergence of Markov chains. There are two components, a theoretical one that utilize an auxiliary function (a solution of a certain equation) to bound convergence, and an empirical one that utilize deep neural networks to approximate the auxiliary function. Furthermore, the authors provide statistical guarantees on the approximation.
strength_weakness: strength_weakness: I'm not able to find effective weaknesses in this paper.
questions: questions: No.
limitations: limitations: I'm not seeing any limitations of essence not addressed in this paper.
paper_summary: paper_summary: The authors proposed a novel computational method to estimate the convergence rate of general Markov Chains. They utilized neural network to verify if contration drift (CD) holds for a given Markov Chain. As an extension to a prior work [1], the authors provided further theoretical analysis of their methods, and proposed an explicit formula of the convergence rate. Various numerical experiments were conducted to justify the applicability.
strength_weakness: strength_weakness: 1. The paper is somewhat intricate and difficult follow. I notice that it is built upon [1], but more detailed background and related work should be provided to improve readability. Specifically, the authors should discuss more previous methods if there exist and make comparisons theoretically or numemrically. 

2. In experiments, the authors only showed the results of their algorithms. I think the convergence rate should be verifed through directly simulating the Markov chain and computing Wasserstein distance between the sampled distribution and stationary distribution.

3. The authors only conducted experiments in 2D, lacking high dimensional examples.

4. A more fundamental problem is that how to verify whether CD holds in practice in high dimensional case. The inequality is defined in the pointwise sense, which is difficult to verify due to the curse of dimensionality.

[1] Qu, Yanlin, Jose Blanchet, and Peter Glynn. "Computable Bounds on Convergence of Markov Chains in Wasserstein Distance." arXiv preprint arXiv:2308.10341 (2023).
questions: questions: 1. What is the relationship between CD and conventional inequailities such as Poincare inequality, Log Sobolev inequality?

1. I wonder the choices of $U$. In all the experiments, $U$ is fixed as a constant. How is this constant selected? Also, is there any other choice of $U$ that can be considered in practice? How does different choices of $U$ lead to different convergence rate estimation? What is the practicability of sequentially learning neural network in section 3.4?

2. The convergence rate in Theorem 3 depends on $\inf/\sup V$. How is inf/sup $V_\theta$ computed in practice?

3. Does there always exist $U,V$ to satisfy CD? If not, how can one know that from DCDC? In another word, how can one tell the failure of training of network from the fact that CD doesn't hold?
limitations: limitations: The limitations are pointed out in the weakness and question part.

Comments: title: Please read the authors’ rebuttal and reply by August 13, 2024 (11:59 PM, AOE time)
Thank you
comment: Dear Reviewers,

Thank you for your hard work during the review process. The authors have responded to your initial reviews. **If you haven’t already done so, please take the time to carefully read their responses.** It is crucial for the authors to receive acknowledgment of your review by the deadline for the author-reviewer discussion period, which is August 13, 2024 (11:59 PM, Anywhere on Earth). Please address any points of disagreement with the authors as early as possible.

Best,

Your AC
Thank you for your reply. The suggestion of trying easier Markov chains for exact rates is very helpful and we will include this comparison. Since we mainly care about whether our bounds are good for reasonably large $n$’s (but not the first several $n$’s), given the exponential convergence, we may need to estimate some very small $W(X_n,X_\infty)$, which is pretty hard in general, but we could do it efficiently in a non-trivial multidimensional process, using importance sampling.
Thanks for you feedback addressing most of my concerns. I still believe some more experiments mentioned in W2 can be conducted to strengthen the paper. I understand the expensive cost of simulation a general non-trivial Markov chain, but you can try some easier one. For instance, in the SGD for logistic regression case, the computational cost is low and thus one can run sufficienly many steps (e.g. 100k) to get sufficiently accurate estimation of $X_\infty$. Then the true convergence rate can be computed (at least approximated) and compared with the theoretical results. Nevertheless, I think this work is indeed interesting and constructive after reading the author rebuttal. I'm happy to raise my score.
Thank you for your positive feedback on my positive feedback. I'll read the other rebuttals.