Explicit Mean-Square Error Bounds
for Monte-Carlo and Linear Stochastic Approximation

As a prototypical example of stochastic approximation, Markov chain Monte Carlo (MCMC) proceeds by constructing an ergodic Markov chain $\Phi$ with invariant measure $\pi$ so as to estimate $\pi(F)=\int F(z)\,\pi(dz)$ for some function $F:{\sf Z}\rightarrow\mathbb{R}^{d}$. One then simulates $\Phi_{1},\Phi_{2},\dots,\Phi_{n+1}$ to obtain the estimates

$$\theta_{n+1}=\frac{1}{n+1}\sum_{k=1}^{n+1}F(\Phi_{k})$$

(8)

This is an instance of the SA recursion (1):

$$\theta_{n+1}=\theta_{n}+\frac{1}{n+1}(-\theta_{n}+F(\Phi_{n+1}))$$

(9)

Subtracting $\theta^{*}=\pi(F)$ from both sides of (9) gives, with $\widetilde{\theta}_{n}=\theta_{n}-\pi(F)$,

$$\widetilde{\theta}_{n+1}=\widetilde{\theta}_{n}+\frac{1}{n+1}(-\widetilde{\theta}_{n}+F(\Phi_{n+1})-\pi(F))$$

which is (3) in a special case: $A=-I$, $\Delta_{n+1}=F(\Phi_{n+1})-\pi(F)$ and $\alpha_{n}=1/n$.

A significant part of the literature on MCMC focuses on finding Markov chains whose marginals approach the invariant measure $\pi$ quickly. Error estimates for MCMC have only been studied under rather restrictive settings. For instance, under the assumption of uniform ergodicity of $\Phi$ and uniform boundedness of $F$ (which rarely hold in practice outside of a finite state space), a generalized Hoeffding’s inequality was obtained in [19] to obtain the PAC-style error bound (7). We can not expect Hoeffding’s bound if either of these assumptions is relaxed. Consider the simplest countable state space Markov chain: the M/M/1 queue with uniformization, defined with ${\sf Z}=\{0,1,2,\dots\}$ and

$$\Phi_{n+1}=\begin{cases}\Phi_{n}+1&\text{prob.\ $\alpha$}\\ \max(\Phi_{n}-1,0)&\text{prob.\ $\mu=1-\alpha$}\end{cases}$$

This is a reversible, geometrically ergodic Markov chain when $\rho=\alpha/\mu<1$, with geometric invariant distribution $\pi(z)=(1-\rho)\rho^{z}$, $z\geq 0$. It is shown in [28] that the error bound (7) fails for most unbounded functions $F$. The question is looked at in greater depth in [17, 16], where asymptotic bounds are obtained for the special case $F(z)\equiv z$. An asymptotic version of (7) is obtained for the lower tail:

$\displaystyle\lim_{n\to\infty}\frac{1}{n}\log\Bigl{(}{\sf P}\{\widetilde{\theta}_{n}\leq-\varepsilon\}\Bigr{)}$

$\displaystyle=-I_{0}(\varepsilon)$

(10)

in which the right hand side is strictly negative and finite valued for positive $\varepsilon$ in a neighborhood of zero. An entirely different scaling is required for the upper tail:

$\displaystyle\lim_{n\to\infty}\frac{1}{n}\log\Bigl{(}{\sf P}\Bigl{\{}\frac{\widetilde{\theta}_{n}}{n}\geq\varepsilon\Bigr{\}}\Bigr{)}$

$\displaystyle=-J_{0}(\varepsilon)$

(11)

where again the right hand side is strictly negative and finite valued for $\varepsilon>0$ sufficiently small. It follows from (11) that the PAC-style bound (7) is not attainable.

The theory of this paper also applies to TD-learning. In this case, the Markov chain $\Phi$ contains as one component a state process for a system to be controlled.

Consider a Markov chain $X$ evolving on a (Polish) state space ${\sf X}$. Given a cost function $c:{\sf X}\to\mathbb{R}$, and a discount factor $\beta\in(0,1)$, the goal in TD-learning is to approximate the solution $h:{\sf X}\to\mathbb{R}$ to the Bellman equation:

$$h(x)=c(x)+\beta{\sf E}[h(X_{n+1})\mid X_{n}=x]$$

(12)

This functional equation can be recast:

	$\displaystyle{\sf E}[$	$\displaystyle\mathcal{D}(h,\Phi_{n+1})\mid\Phi_{0}\,\dots\,\Phi_{n}]=0$		(13a)
	where	$\displaystyle\mathcal{D}(h,\Phi_{n+1})\mathbin{:=}c(X_{n})+\beta h(X_{n+1})-h(X_{n})\,,\quad\Phi_{n+1}\mathbin{:=}(X_{n+1},X_{n})\,,\quad n\geq 0$		(13b)

Equation (13a) may be regarded as motivation for the TD-learning algorithms of [36, 38].

Consider a linearly parameterized family of candidate approximations $\{h^{\theta}(x)=\theta^{\intercal}\psi(x):\theta\in\mathbb{R}^{d}\}$, where $\psi:{\sf X}\to\mathbb{R}^{d}$ denotes the $d$ basis functions. The goal in TD-learning is to solve the Galerkin relaxation of (13a,13b):

$${\sf E}[\mathcal{D}(h^{\theta^{*}},\Phi_{n+1})\zeta_{n}]=0$$

(14)

where $\{\zeta_{n}\}$ is a $d$-dimensional stochastic process, adapted to $\Phi$, and the expectation is with respect to the steady state distribution. In particular, $\zeta_{n}\equiv\psi(X_{n})$ in TD($0$) learning, so that the goal in this case is to find $\theta^{*}\in\mathbb{R}^{d}$ such that:

	$\displaystyle{\sf E}$	$\displaystyle\big{[}\mathcal{D}(h^{\theta^{*}},\Phi_{n+1})\psi(X_{n})\big{]}=0$		(15)
		$\displaystyle\mathcal{D}(h^{\theta^{*}},\Phi_{n+1})=c(X_{n})+\beta h^{\theta^{*}}(X_{n+1})-h^{\theta^{*}}(X_{n})$

The TD($0$) algorithm is the SA recursion (1) applied to solve (15):

	$\displaystyle\theta_{n+1}$	$\displaystyle=\theta_{n}+\alpha_{n+1}d_{n+1}\psi(X_{n})$		(16)
	$\displaystyle d_{n+1}$	$\displaystyle=c(X_{n})+\beta h^{\theta_{n}}(X_{n+1})-h^{\theta_{n}}(X_{n})$

Denoting

	$\displaystyle A_{n+1}$	$\displaystyle\mathbin{:=}\psi(X_{n})\big{(}\beta\psi(X_{n+1})-\psi(X_{n})\big{)}^{\intercal}$
	$\displaystyle b_{n+1}$	$\displaystyle\mathbin{:=}-c(X_{n})\psi(X_{n})$

the algorithm (16) can be rewritten as:

$$\theta_{n+1}=\theta_{n}+\alpha_{n+1}\big{(}A_{n+1}\theta_{n}-b_{n+1}\big{)}$$

(17)

Note that $\theta^{*}$ from (14) solves the linear equation ${\sf E}[A_{n+1}]\theta^{*}={\sf E}[b_{n+1}]$. Subtracting $\theta^{*}$ from both sides of (17) gives, with $\widetilde{\theta}_{n}=\theta_{n}-\theta^{*}$,

$$\widetilde{\theta}_{n+1}=\widetilde{\theta}_{n}+\alpha_{n+1}[A\widetilde{\theta}_{n}+(A_{n+1}-A)\widetilde{\theta}_{n}+A_{n+1}\theta^{*}-b_{n+1}]$$

(18)

Under mild conditions, we show through coupling that iteration (18) can be closely approximated by the linear SA recursion (3) with matrix $A={\sf E}[A_{n+1}]$ and noise sequence $\Delta_{n+1}=A_{n+1}\theta^{*}-b_{n+1}$. In particular, the two recursions have the same asymptotic covariance if the matrix ${\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}I+A$ is Hurwitz (see Section 2.3).

Under general assumptions on the Markov chain $X$, and the basis functions $\psi$, it is known that matrix $A=\partial{\sf E}[\mathcal{D}(h^{\theta^{*}},\Phi_{n+1})]={\sf E}[A_{n+1}]$ is Hurwitz, and that the sequence of estimates $\{\theta_{n}\}$ converges to $\theta^{*}$ [38]. However, when the discount factor $\beta$ is close to $1$, it can be shown that $\lambda_{\max}>-{\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}$ (where $\lambda_{\max}$ denotes the largest eigenvalue of $A$), and is in fact close to $0$ under mild additional assumptions [14, 11, 13]. It follows that the algorithm has infinite asymptotic covariance: full details and finer results can be deduced from Theorems 2.4 and 2.7.

The SNR algorithm is defined as follows:

	$\displaystyle\theta_{n+1}$	$\displaystyle=\theta_{n}-\alpha_{n+1}\widehat{A}_{n+1}^{-1}d_{n+1}\psi(X_{n})$		(19)
	$\displaystyle d_{n+1}$	$\displaystyle=c(X_{n})+\beta h^{\theta_{n}}(X_{n+1})-h^{\theta_{n}}(X_{n})$
	$\displaystyle\widehat{A}_{n+1}$	$\displaystyle=\widehat{A}_{n}+\alpha_{n+1}\big{[}A_{n+1}-\widehat{A}_{n}\big{]}$		(20)

Under the assumption that the sequence of matrices $\{\widehat{A}_{n}:n\geq 0\}$ is invertible for each $n$, it is shown in [14, 11] that the sequence of estimates obtained using (19,20) are identical to the parameter estimates obtained using the LSTD($0$) algorithm: $\theta_{n+1}=\widehat{A}_{n+1}^{-1}\widehat{b}_{n+1}$, with

	$\displaystyle\widehat{A}_{n+1}$	$\displaystyle=\widehat{A}_{n}+\alpha_{n+1}\big{[}A_{n+1}-\widehat{A}_{n}\big{]}$
	$\displaystyle\widehat{b}_{n+1}$	$\displaystyle=\widehat{b}_{n}+\alpha_{n+1}\big{[}b_{n+1}-\widehat{b}_{n}\big{]}$

Consequently, the LSTD($0$) algorithm achieves the optimal asymptotic covariance.

Q-learning and many other RL algorithms can also be cast as SA recursions. They are no longer linear, but it is anticipated that bounds can be obtain in future research through linearization [18].

Finite time performance bounds for linear stochastic approximation were obtained in many prior papers, subject to the assumption that the noise sequence $\{\Delta_{n}\}$ appearing in (3) is a martingale difference sequence [10, 26]. This assumption is rarely satisfied in the applications of interest to the authors.

Much of the literature on finite time bounds for linear SA recursions with Markovian noise has been recent. For constant step-size algorithms with step-size $\alpha$, it follows from analysis in [6] that the pair process $(\theta_{n},\Phi_{n})$ is a geometrically ergodic Markov chain, and the covariance of $\theta_{n}$ is $O(\alpha)$ in steady state. Finite time bounds of order $O(\alpha)$ were obtained in [37, 3, 35, 21]. Unfortunately, these bounds are not tight, and hence their value for algorithm design is limited.

Mean-square error bounds have also been obtained for diminishing step-size algorithms, to establish the optimal rate of convergence ${\sf E}[\|\widetilde{\theta}_{n}\|^{2}]\leq b_{\theta}/n$ [35, 3, 8]. The constant $b_{\theta}$ is a function of the mixing time of the underlying Markov chain. These results require strong assumptions (uniform ergodicity of the Markov chain), and do not obtain the optimal constant $b^{*}_{\theta}=\hbox{\rm trace\,}(\Sigma_{\theta})$. Rather than parameter estimation error, finite time bounds are obtained in [22] for ${\sf E}[\|\bar{f}(\theta_{n})\|^{2}]$, which may be regarded as a far more relevant performance criterion. Bounds are obtained for Markovian models, subject to the existence of a Lyapunov function similar to what is assumed in the present work. It is again not clear if the resulting bounds are tight, or have value in algorithm design.

The main contribution of this paper is a general framework for analyzing the finite time performance of linear stochastic approximation algorithms with Markovian noise, and vanishing step-size (required to achieve the optimal rate of convergence of Chung-Ruppert-Polyak). The M/M/1 queue example illustrates plainly that Markovian noise introduces challenges not seen in the “white noise” setting, and that the finite-$n$ error bound (7) cannot be obtained without substantial restrictions. Even under the assumptions of [19] (uniform ergodicity, and bounded noise), the resulting bounds are extremely loose and hence may give little insight for algorithm design. Our approach allows us to obtain explicit bounds, without the uniform boundedness assumption of noise that is frequently imposed in the literature [3, 35, 21, 8]. Instead, it is assumed that the Markovian noise is $V$-uniform ergodic; an assumption that is far weaker than geometric or uniform mixing.

Our starting point is the classical martingale approximation of the noise used in CLT analysis of Markov chains [29, Chapter 17] and used in the analysis of SA recursions since Metivier and Priouret [27]. Under mild assumptions on the Markov chain, each $\Delta_{n}$ can be expressed as the sum of a martinagle difference and a telescoping term. The solution of the linear recursion (3) is decomposed as a sum of the respective responses:

$$\widetilde{\theta}_{n}=\widetilde{\theta}_{n}^{\mathcal{M}}+\widetilde{\theta}_{n}^{\mathcal{T}}$$

(21)

The challenge is to obtain explicit bounds on the mean square error for each term.

We say that a deterministic vector-valued sequence $\{e_{n}\}$ converges to zero at rate $1/n^{\varrho_{0}}$ if

$$\lim_{n\to\infty}n^{\varrho}\|e_{n}\|=\begin{cases}0,&\text{if }\;\varrho<\varrho_{0}\\ \infty,&\text{if }\;\varrho>\varrho_{0}\end{cases}$$

Bounds for the mean-square error are obtained in Thm. 2.4, subject to conditions on both the matrix $A$ and the noise sequence. In summary, under general assumptions on $\{\Delta_{n}\}$,

• (i)

  The bound (5) holds if ${\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}I+A$ is Hurwitz.

• (ii)

  If $I+A$ is Hurwitz, then the finer bound (6) holds.

• (iii)

  If there is an eigenvalue of $A$ satisfying $\text{Real}(\lambda)>-{\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}$, and corresponding left-eigenvalue $v$ that lies outside of the null-space of the asymptotic covariance of the noise sequence, then

  $$\lim_{n\to\infty}n^{2\rho}{\sf E}[|v^{\intercal}\widetilde{\theta}_{n}|^{2}]=\begin{cases}0,&\quad\varrho<\varrho_{0}\\ \infty,&\quad\varrho>\varrho_{0}\end{cases}$$

  (22)

  with $\rho_{0}=|\text{Real}(\lambda)|$. The convergence of ${\sf E}[\|\widetilde{\theta}_{n}\|^{2}]$ to zero is thus no faster than $n^{-2\rho_{0}}$.

Consider the linear SA recursion (3), with the noise sequence $\{\Delta_{n}\}$ defined in (2). We use the following notation to explicitly represent the noise as a function of $\Phi_{n}$:

$$f^{*}(\Phi_{n})\mathbin{:=}\Delta_{n}=f(\theta^{*}\,,\Phi_{n})$$

(23)

A form of geometric ergodicity is assumed throughout. To apply standard theory, it is assumed that the state space ${\sf Z}$ is Polish (the standing assumption in [29]). We fix a measurable function $V\colon{\sf Z}\to[1,\infty)$, and let $L_{\infty}^{V}$ denote the set of measurable functions $g\colon{\sf Z}\to\mathbb{R}$ satisfying

$$\|g\|_{V}\mathbin{:=}\sup_{z\in{\sf Z}}\frac{|g(z)|}{V(z)}<\infty$$

The Markov chain $\Phi$ is assumed to be $V$-uniformly ergodic: there exists $\rho\in(0,1)$, and $B_{V}<\infty$ such that for each $g\in L_{\infty}^{V}$, $z\in{\sf Z}$,

$$\Big{|}{\sf E}[g(\Phi_{n})\mid\Phi_{0}=z]-\pi(g)\Big{|}\leq B_{V}\|g\|_{V}\rho^{n}V(z)\,,\quad n\geq 0$$

(24)

where $\pi$ is the unique invariant measure, and $\pi(g)=\int g(z)\,\pi(dz)$ is the steady state mean.

The uniform bound (24) is not a strong assumption. For example, it is satisfied for the M/M/1 queue described in Section 1.1 with $V(z)=\exp(\varepsilon_{0}z)$, for $\varepsilon_{0}>0$ sufficiently small, with $z\in{\sf Z}=\{0,1,\dots\}$ [29, Thm. 16.4.1].

The following are imposed throughout:

We now explain the decomposition (21). Each of the two sequences $\{\widetilde{\theta}_{n}^{\mathcal{M}}\}$ and $\{\widetilde{\theta}_{n}^{\mathcal{T}}\}$ evolves as a stochastic approximation sequence, differentiated by the inputs and initial conditions:

	$\displaystyle\widetilde{\theta}_{n+1}^{\mathcal{M}}$	$\displaystyle=\widetilde{\theta}_{n}^{\mathcal{M}}+\alpha_{n+1}\bigl{[}A\widetilde{\theta}_{n}^{\mathcal{M}}\!+\Delta_{n+2}^{m}\bigr{]}\,,$	$\displaystyle\widetilde{\theta}_{0}^{\mathcal{M}}=\widetilde{\theta}_{0}$		(30a)
	$\displaystyle\widetilde{\theta}_{n+1}^{\mathcal{T}}$	$\displaystyle=\widetilde{\theta}_{n}^{\mathcal{T}}+\alpha_{n+1}\bigl{[}A\widetilde{\theta}_{n}^{\mathcal{T}}\!+\!Z_{n+1}-Z_{n+2}\bigr{]}\,,$	$\displaystyle\widetilde{\theta}_{0}^{\mathcal{T}}=0$		(30b)

The second recursion admits a more tractable realization through the change of variables, $\Xi_{n}=\widetilde{\theta}_{n}^{\mathcal{T}}+\alpha_{n}Z_{n+1}$, $n\geq 1$.

Lemma 2.2 combined with (30) gives

$$\widetilde{\theta}_{n}=\widetilde{\theta}_{n}^{(1)}+\widetilde{\theta}_{n}^{(2)}+\widetilde{\theta}_{n}^{(3)}$$

(32)

where $\widetilde{\theta}_{n}^{(1)}\!=\widetilde{\theta}_{n}^{\mathcal{M}}$, $\widetilde{\theta}_{n}^{(2)}\!=\Xi_{n}$, and $\widetilde{\theta}_{n}^{(3)}\!=\!-\alpha_{n}Z_{n+1}$ for $n\geq 1$. Note that $\widetilde{\theta}_{n}^{\mathcal{T}}=\widetilde{\theta}_{n}^{(2)}+\widetilde{\theta}_{n}^{(3)}$.

It is more convenient to work directly with the recursion for the scaled sequence:

Denote $\widetilde{\theta}_{n}^{\varrho,(i)}=n^{\varrho}\widetilde{\theta}_{n}^{(i)}$ for each $i$. Lemma 2.3 combined with (32) gives

$$\widetilde{\theta}^{\varrho}_{n}=\widetilde{\theta}_{n}^{\varrho,(1)}+\widetilde{\theta}_{n}^{\varrho,(2)}+\widetilde{\theta}_{n}^{\varrho,(3)}$$

(34)

The first two sequences evolve as SA recursions:

	$\displaystyle\widetilde{\theta}_{n+1}^{\varrho,(1)}=\widetilde{\theta}_{n}^{\varrho,(1)}+$	$\displaystyle\alpha_{n+1}\bigl{[}[\varrho_{n}I+A(n,\varrho)]\widetilde{\theta}_{n}^{\varrho,(1)}+(n+1)^{\varrho}\Delta_{n+2}^{m}\bigr{]}\,,$	$\displaystyle\widetilde{\theta}_{0}^{\varrho,(1)}\!=\!\widetilde{\theta}^{\varrho}_{0}$		(35a)
	$\displaystyle\widetilde{\theta}_{n+1}^{\varrho,(2)}=\widetilde{\theta}_{n}^{\varrho,(2)}+$	$\displaystyle\alpha_{n+1}\bigl{[}[\varrho_{n}I+A(n,\varrho)]\widetilde{\theta}_{n}^{\varrho,(2)}-(n+1)^{\varrho}\alpha_{n}[I+A]Z_{n+1}\bigr{]}\,,$	$\displaystyle\widetilde{\theta}^{\varrho,(2)}_{1}\!=\!\Xi_{1}$		(35b)

Fix the initial condition $(\Phi_{0},\widetilde{\theta}_{0})$, and denote $\text{\rm Cov}\,(\theta_{n})={\sf E}[\widetilde{\theta}_{n}\widetilde{\theta}_{n}^{\intercal}]$ and $\Sigma_{Z}={\sf E}_{\pi}[Z_{n}Z_{n}^{\intercal}]$. The following summarizes bounds on the convergence rate of ${\sf E}[\|\widetilde{\theta}_{n}\|^{2}]=\hbox{\rm trace\,}(\text{\rm Cov}\,(\theta_{n}))$.

The proof of Thm. 2.4 is contained in Section 2.5. The following negative result is a direct corollary of Thm. 2.4 (ii):

One challenge in extension to nonlinear recursions is that the noise sequence depends on the parameter estimates (recall (2)). This is true even for TD learning with linear function approximation (see (18) and surrounding discussion). Extension to these recursions is obtained through coupling.

Consider the error sequence for a random linear recursion

$$\widetilde{\theta}^{\circ}_{n+1}=\widetilde{\theta}^{\circ}_{n}+\alpha_{n+1}[A_{n+1}\widetilde{\theta}^{\circ}_{n}+A_{n+1}\theta^{*}-b_{n+1}]$$

(36)

subject to the following assumptions:

• (A4)

  The sequences $\{A_{n},b_{n}\}$ are functions of the Markov chain:

  $$A_{n}=\mathcal{A}(\Phi_{n})\,,\quad b_{n}=\mathcal{B}(\Phi_{n})\,,$$

  which satisfy $\|\mathcal{A}_{i,j}^{2}\|_{V}<\infty$, $\|\mathcal{B}_{i}^{2}\|_{V}<\infty$ for each $1\leq i,j\leq d$. The steady state means are denoted $A={\sf E}_{\pi}[A_{n}]$, $b={\sf E}_{\pi}[b_{n}]$. Moreover, the matrix $A$ is Hurwitz, and $\theta^{*}=A^{-1}b$.

The proof of the theorem is via coupling with (3). For this we write (36) in the suggestive form

$$\widetilde{\theta}^{\circ}_{n+1}=\widetilde{\theta}^{\circ}_{n}+\alpha_{n+1}[A\widetilde{\theta}^{\circ}_{n}+\Delta_{n+1}+(A_{n+1}-A)\widetilde{\theta}^{\circ}_{n}]\,,\qquad\Delta_{n+1}=A_{n+1}\theta^{*}-b_{n+1}$$

(37)

With common initial condition $\Phi_{0}$, the sequence $\{\widetilde{\theta}^{\circ}_{n}\}$ is compared with the error sequence $\{\widetilde{\theta}^{\bullet}_{n}\}$ for the corresponding linear SA algorithm:

$$\widetilde{\theta}^{\bullet}_{n+1}=\widetilde{\theta}^{\bullet}_{n}+\alpha_{n+1}[A\widetilde{\theta}^{\bullet}_{n}+\Delta_{n+1}]$$

The difference sequence $\{\mathcal{E}_{n}\mathbin{:=}\widetilde{\theta}^{\circ}_{n}-\widetilde{\theta}^{\bullet}_{n}\}$ evolves according to (3), but with a vanishing noise sequence:

$$\mathcal{E}_{n+1}=\mathcal{E}_{n}+\alpha_{n+1}[A\mathcal{E}_{n}+(A_{n+1}-A)\widetilde{\theta}^{\circ}_{n}]$$

(38)

By decomposing $A_{n+1}-A$ into martingale difference and telescoping sequences based on Poisson’s equation, the technique used to prove Thm. 2.4 can be used to obtain the following bound on the mean-square coupling error.

Let $\lambda=-\varrho_{0}+ui$ denote an eigenvalue of the matrix $A$ with largest real part (i.e., $\varrho_{0}$ is minimal).

Thm. 2.7 provides a remarkable bound when $\rho_{0}>1$: it immediately implies Thm. 2.6 because the mean-square coupling error ${\sf E}[\|\mathcal{E}_{n}]\|^{2}$ tends to zero at rate no slower than $n^{-2}$, which is far faster than ${\sf E}[\|\widetilde{\theta}^{\bullet}_{n}\|^{2}]\approx\hbox{\rm trace\,}(\Sigma_{\theta})n^{-1}$ (implied by Thm. 2.4).

An alert reader may observe that Theorems 2.6 and 2.7 leave out a special case: consider $\rho_{0}<{\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}$, so that the rate of convergence of ${\sf E}[\|\widetilde{\theta}^{\bullet}_{n}\|^{2}]$ is the sub-optimal value $n^{-2\rho_{0}}$. The bound obtained in Thm. 2.7 remains valuable, in the sense that it combined with Thm. 2.4 (ii) implies the rate of convergence of ${\sf E}[\|\widetilde{\theta}^{\circ}_{n}]\|^{2}$ is no slower than $n^{-2\rho_{0}}$. However, because ${\sf E}[\|\mathcal{E}_{n}]\|^{2}$ and ${\sf E}[\|\widetilde{\theta}^{\bullet}_{n}\|^{2}]$ tend to zero at the same rate, we cannot rule out the possibility that $\widetilde{\theta}^{\circ}_{n}=\mathcal{E}_{n}+\widetilde{\theta}^{\bullet}_{n}$ converges to zero much faster. In particular, it remains to prove that if there is an eigenvalue $\lambda$ of $A$ that satisfies $-\varrho_{0}=\text{Real}(\lambda)>-{\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}$, and an eigenvector $v\neq 0$ satisfying $\Sigma_{\Delta}v\neq 0$, then, ${\sf E}[|v^{\intercal}\widetilde{\theta}_{n}^{\circ}|^{2}]$ converges to $0$ at rate $n^{-2\varrho_{0}}$.

Thm. 2.4 indicates that the convergence rate of $\text{\rm Cov}\,(\theta_{n})$ is determined jointly by the matrix $A$, and the martingale difference component of the noise sequence $\{\Delta_{n}\}$. Convergence of $\{\widetilde{\theta}_{n}\}$ can be slow if the matrix $A$ has eigenvalues close to zero.

The result also explains the slow convergence of some reinforcement learning algorithms. For instance, the matrix $A$ in Watkins’ Q-learning has at least one eigenvalue with real part greater than or equal to $-(1-\beta)$, where $\beta$ is the discount factor appearing in the Markov decision process [40, 14, 11]. Since $\beta$ is usually close to one, Thm. 2.4 implies that the convergence rate of the algorithm is much slower than $n^{-1}$. Under the assumption that $A$ is Hurwitz, the $1/n$ convergence rate is guaranteed by the use of a modified step-size sequence $\alpha_{n}=g/n$, with $g>0$ chosen so that the matrix ${\mathchoice{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{1}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}{\genfrac{}{}{}{3}{1}{2}}}I+gA$ is Hurwitz. Corollary 2.8 follows directly from Thm. 2.4 (i).

We can also ensure the $1/n$ convergence rate by using a matrix gain. Provided $A$ is invertible, and if it is known beforehand, $\alpha_{n}=-A^{-1}/n$ is the optimal step-size sequence (in terms of minimizing the asymptotic covariance) [1, 25, 13]. The SQNR algorithm of [33] and the Zap-SNR algorithm [14, 11] provide general approaches to recursively estimate the optimal matrix gain.

The next subsection is dedicated to the proof of Thm. 2.4. The proofs of the technical results are contained in the Appendix A.

Denote $\text{\rm Cov}\,(\theta_{n}^{(i)})={\sf E}[\widetilde{\theta}_{n}^{(i)}(\widetilde{\theta}_{n}^{(i)})^{\intercal}]$ and $\Sigma_{n}^{\varrho,(i)}={\sf E}[\widetilde{\theta}^{\varrho,(i)}(\widetilde{\theta}^{\varrho,(i)})^{\intercal}]=n^{2\varrho}\text{\rm Cov}\,(\theta_{n}^{(i)})$ for each $i$ in (34). The proof proceeds by establishing the convergence rate for each $\text{\rm Cov}\,(\theta_{n}^{(i)})$. The main challenges are the first two: $\text{\rm Cov}\,(\theta_{n}^{(1)})$ and $\text{\rm Cov}\,(\theta_{n}^{(2)})$, for which explicit bounds are obtained by studying recursions of the scaled sequences. Bounding $\widetilde{\theta}_{n}^{{(3)}}=-\alpha_{n}Z_{n+1}$ is trivial.