Non-asymptotic convergence analysis for the Unadjusted Langevin Algorithm
: In this paper, we study a method to sample from a target distribution over having a positive density with respect to the Lebesgue measure, known up to a normalisation factor. This method is based on the Euler discretization of the overdamped Langevin stochastic differential equation associated with . For both constant and decreasing step sizes in the Euler discretization, we obtain non-asymptotic bounds for the convergence to the target distribution in total variation distance. A particular attention is paid to the dependency on the dimension , to demonstrate the applicability of this method in the high dimensional setting. These bounds improve and extend the results of [dalalyan:2014].
lemmatheorem \aliascntresetthelemma \newaliascntcorollarytheorem \aliascntresetthecorollary \newaliascntpropositiontheorem \aliascntresettheproposition \newaliascntdefinitiontheorem \aliascntresetthedefinition \newaliascntremarktheorem \aliascntresettheremark
Non-asymptotic convergence analysis for the ULA
class=AMS] \kwd[primary ]65C05, 60F05, 62L10 \kwd[; secondary ]65C40, 60J05,93E35
total variation distance\kwdLangevin diffusion \kwdMarkov Chain Monte Carlo \kwdMetropolis Adjusted Langevin Algorithm \kwdRate of convergence
Sampling distributions over high-dimensional state-spaces is a problem which has recently attracted a lot of research efforts in computational statistics and machine learning (see [cotter:roberts:stuart:white:2013] and [andrieu:defreitas:doucet:jordan:2003] for details); applications include Bayesian non-parametrics, Bayesian inverse problems and aggregation of estimators. All these problems boil down to sample a target distribution having a density w.r.t. the Lebesgue measure on , known up to a normalisation factor where is continuously differentiable. We consider a sampling method based on the Euler discretization of the overdamped Langevin stochastic differential equation (SDE)
where is a -dimensional Brownian motion. It is well-known that the Markov semi-group associated with the Langevin diffusion is reversible w.r.t. . Under suitable conditions, the convergence to takes place at geometric rate. Precise quantitative estimates of the rate of convergence with explicit dependency on the dimension of the state space have been recently obtained using either functional inequalities such as Poincaré and log-Sobolev inequalities (see [bakry:cattiaux:guillin:2008, cattiaux:guillin:2009] [bakry:gentil:ledoux:2014]) or by coupling techniques (see [eberle:2015]). The Euler-Maruyama discretization scheme associated to the Langevin diffusion yields the discrete time-Markov chain given by
where is an i.i.d. sequence of standard Gaussian -dimensional random vectors and is a sequence of step sizes, which can either be held constant or be chosen to decrease to . The idea of using the Markov chain to sample approximately from the target has been first introduced in the physics literature by [parisi:1981] and popularised in the computational statistics community by [grenander:1983] and [grenander:miller:1994]. It has been studied in depth by [roberts:tweedie:1996], which proposed to use a Metropolis-Hastings step at each iteration to enforce reversibility w.r.t. leading to the Metropolis Adjusted Langevin Algorithm (MALA). They coin the term unadjusted Langevin algorithm (ULA) when the Metropolis-Hastings step is skipped.
The purpose of this paper is to study the convergence of the ULA algorithm. The emphasis is put on non-asymptotic computable bounds; we pay a particular attention to the way these bounds scale with the dimension and constants characterizing the smoothness and curvature of the potential . Our study covers both constant and decreasing step sizes and we analyse both the ”finite horizon” (where the total number of simulations is specified before running the algorithm) and ”any-time” settings (where the algorithm can be stopped after any iteration).
When the step size is constant, under appropriate conditions (see [roberts:tweedie:1996]), the Markov chain is -uniformly geometrically ergodic with a stationary distribution . With few exceptions, the stationary distribution is different from the target . If the step size is small enough, then the stationary distribution of this chain is in some sense close to . We provide non-asymptotic bounds of the -total variation distance between and , with explicit dependence on the step size and the dimension . Our results complete and extend the recent works by [dalalyan:tsybakov:2012] and [dalalyan:2014].
When decreases to zero, then is a non-homogeneous Markov chain. If in addition , we show that the marginal distribution of this non-homogeneous chain converges, under some mild additional conditions, to the target distribution , and provide explicit bounds for the convergence. Compared to the related works by [lamberton:pages:2002], [lamberton:pages:2003], [lemaire:2005] and [lemaire:menozzi:2010], we establish not only the weak convergence of the weighted empirical measure of the path to the target distribution but a much stronger convergence in total variation, similarly to [dalalyan:2014], where the strongly log-concave case is considered.
The paper is organized as follows. In Section 2, the main convergence results are stated under abstract assumptions. We then specialize in Section 3 these results to different classes of densities. The proofs are gathered in Section 4. Some general convergence results for diffusions based on reflection coupling, which are of independent interest, are stated in Section 5.
Notations and conventions
denotes the Borel -field of and the set of all Borel measurable functions on . For set . Denote by the space of finite signed measure on and . For and a -integrable function, denote by the integral of w.r.t. . Let be a measurable function. For , the -norm of is given by . For , the -total variation distance of is defined as
If , then is the total variation denoted by .
For , denote by the set of measurable functions such that . For , the variance of under is denoted by . For all functions such that , the entropy of with respect to is defined by
Let and be two probability measures on . If , we denote by the Radon-Nikodym derivative of w.r.t. . Denote for all by the scalar product of and and the Euclidean norm of . For , denote by , the set of -times continuously differentiable functions . For , denote by the gradient of and the Laplacian of . For all and , we denote by , the ball centered at of radius . Denote for , the oscillation of a function in the ball by . Denote the oscillation of a bounded function on by . In the sequel, we take the convention that and , for , .
2 General conditions for the convergence of ULA
In this section, we derive a bound on the convergence of the ULA to the target distribution when the Langevin diffusion is geometrically ergodic and the Markov kernel associated with the EM discretization satisfies a Foster-Lyapunov drift inequality.
Consider the following assumption on the potential :
The function is continuously differentiable on and gradient Lipschitz, i.e. there exists such that for all ,
Under L 1, by [ikeda:watanabe:1989, Theorem 2.4-3.1] for every initial point , there exists a unique strong solution to the Langevin SDE (1). Define for all , and , . The semigroup is reversible w.r.t. , and hence admits as its (unique) invariant distribution. In this section, we consider the case where is geometrically ergodic, i.e. there exists such that for any initial distribution and ,
for some constant . Denote by the generator associated with the semigroup , given for all by
A twice continuously differentiable function is a Lyapunov function for the generator if there exist , and such that,
By [roberts:tweedie:1996, Theorem 2.2], if in (4) is a non-empty compact set, then the Langevin diffusion is geometrically ergodic.
Consider now the EM discretization of the diffusion (2). Let be a sequence of positive and nonincreasing step sizes and for , denote by
For , consider the Markov kernel given for all and by
The discretized Langevin diffusion given in (2) is a time-inhomogeneous Markov chain, for and , where and
with the convention that for , , is the identity operator. Under L 1, the Markov kernel is strongly Feller, irreducible and strongly aperiodic. We will say that a function satisfies a Foster-Lyapunov drift condition for if there exist constants , and such that, for all
The particular form of (6) reflects how the mixing rate of the Markov chain depends upon the step size . If , then for and . A Markov chain with transition kernel is not mixing. Intuitively, as gets larger, then it is expected that the mixing of increases. If for some , satisfies (6), then admits a unique stationary distribution . We use (6) to control quantitatively the moments of the time-inhomogeneous chain. The types of bounds which are needed, are summarised in the following elementary Lemma.
Let . Assume that for all and , (6) holds for some constants and . Let be a sequence of nonincreasing step sizes such that for all . Then for all and , where
The proof is postponed to Section 4.1. ∎
Note that Section 2 implies that where
We give below the main ingredients which are needed to obtain a quantitative bound for for all . This quantity is decomposed as follows: for all ,
To control the first term on the right hand side, we use a method introduced in [dalalyan:tsybakov:2012] and elaborated in [dalalyan:2014]. The second term is bounded using the convergence of the semi-group to , see (3).
The proof follows the same lines as [dalalyan:2014, Lemma 2] but is given for completeness. For , let be the space of continuous functions on taking values in . For all , denote by and the laws on of the Langevin diffusion and of the continuously-interpolated Euler discretization , both started at at time . Denote by the unique strong solution started at at time of the time-inhomogeneous diffusion defined for , by
where for any continuous function and
Girsanov’s Theorem [karatzas:shreve:1991, Theorem 5.1, Corollary 5.16, Chapter 3] shows that and are mutually absolutely continuous and in addition, -almost surely
In the sequel, depending on the conditions on the potential and the techniques of proof, for any given , can have two kinds of upper bounds, either of the form , or , for some function and . In both cases, as shown in Section 2, it is possible to choose as a function of , so that under appropriate conditions on the sequence of step sizes .
There exists such that for all , and . Therefore, we can define for all ,
and . We first show that . The proof goes by contradiction. If we could extract a bounded subsequence . For such sequence, is bounded away from , but which yields to a contradiction. The definition of implies that , showing that
On the other hand, since is nonincreasing, for any ,
The proof follows from (10) using .
For all , define . Note that since , we have . Using and is a nonincreasing sequence, we get for all ,
which shows that the first term in the right side of (10) goes to as goes to infinity. As for the second term, since , we get using that is nonincreasing and ,
Using and , we have , which concludes the proof.
Using (10), we can also assess the convergence of the algorithm for constant step sizes for all . Two different kinds of results can be derived. First, for a given precision , we can try to optimize the step size to minimize the number of iterations required to achieve . Second if the total number of iterations is fixed , we may determine the step size which minimizes .
For , set . Then using the stated expressions of and in (10) concludes the proof. ∎
The proof is a straightforward calculation using (10). ∎
To get quantitative bounds for the total variation distance it is therefore required to get bounds on , and to control . We will consider in the sequel two different approaches to get (3), one based on functional inequalities, the other on coupling techniques. We will consider also increasingly stringent assumptions for the potential . Whereas we will always obtain the same type of exponential bounds, the dependency of the constants on the dimension will be markedly different. In the worst case, the dependency is exponential. It is polynomial when is convex.
3 Practical conditions for geometric ergodicity of the Langevin diffusion and their consequences for ULA
3.1 Superexponential densities
Assume first that the potential is superexponential outside a ball. This is a rather weak assumption (we do not assume convexity here).
The potential is twice continuously differentiable and there exist , and such that for all , , .
The price to pay will be constants which are exponential in the dimension. Under H 1, the potential is unbounded off compact set. Since is continuous, it has a global minimizer , which is a point at which . Without loss of generality, it is assumed that .
The elementary proof is postponed to Section 4.2. ∎
Following [roberts:tweedie:1996, Theorem 2.3], we first establish a drift condition for the diffusion.
The proof, adapted from [roberts:tweedie:1996, Theorem 2.3] and [meyn:tweedie:1993:3, Theorem 6.1], is postponed to Section 4.3. ∎
The proof is postponed to Section 4.4. ∎
The proof is postponed to Section 4.5. ∎
Let , defined for all by . By Section 3.1, is a contraction operator on the space of finite signed measure , , endowed with the norm . It is therefore possible to control . To simplify the notations, we limit our discussion to constant step sizes.
Note that Theorem 2 implies that there exists a constant which does not depend on such that .
It is shown in [talay:tubaro:1991, Theorem 4] that for with polynomial growth, , for some constant , provided that satisfies L 1 and H 1. Our result does not match this bound since . However the bound is uniform over the class of measurable functions satisfying for all , . Obtaining such uniform bounds in total variation is important in Bayesian inference, for example to compute high posterior density credible regions. Our result also strengthens and completes [mattingly:stuart:higham:2002, Corollary 7.5], which states that under H 1 with , for any measurable functions satisfying for all ,
for some , , for some constants and