Bernard Bercu
email: bernard.bercu@math.u-bordeaux1.fr

Pierre Del Moral
email: pierre.del-moral@inria.fr

Arnaud Doucet
email: arnaud@stat.ubc.ca

Fluctuations of Interacting Markov Chain Monte Carlo Models

Keywords: Multivariate and functional central limit theorems, random fields, martingale limit theorems, self-interacting Markov chains, Markov chain Monte Carlo models Théorèmes de la limite centrale multidimensionnels, champs aléatoires, théorèmes limites de martingales, processus en auto-interaction, méthodes de Monte Carlo

We present a functional central limit theorem for a general class of interacting Markov chain Monte Carlo interpretations of discrete generation measure-valued equations. The path space models associated with these stochastic processes belong to the class of nonlinear Markov chains interacting with their empirical occupation measures. We develop an original theoretical analysis based on resolvent operators and semigroup techniques to analyze the fluctuation of their occupation measures around their limiting value. We also present a set of simple regularity conditions that applies to interacting Markov chain Monte Carlo models on path spaces, yielding what seems to be the first fluctuation theorems for this class of self-interacting models.

1 Introduction

Interacting Markov chain Monte Carlo models

Interacting Markov chain Monte Carlo methods (i-MCMC ) are a new class of stochastic simulation algorithms for the numerical solving of a fairly general class of discrete generation nonlinear measure-valued equations. This class of models are defined in the following way. Let (S (l) , S (l) ) l≥0 be a sequence of measurable spaces, and for every l ≥ 0 we denote by P(S (l) ) the set of all probability measures on S (l) . We consider a flow of probability measures π (l) ∈ P(S (l) ) that satisfies an evolution equation of the following form

∀l ≥ 0 Φ l (π (l-1) ) = π (l) (1.1)
for some mappings Φ l : P(S (l-1) ) → P(S (l) ). In this notation system, we simplify the presentation using the convention Φ 0 (π (-1) ) = π (0) , for l = 0. The i-MCMC model associated to (1.1) is defined in terms of a series of Markov transitions. That is, we consider a Markov transition M (0) from S (0) into itself , and a collection of Markov transitions M

µ from S (l) into itself, indexed by the parameter l ≥ 0 and the set of probability measures µ ∈ P(S (l-1) ). We further assume that the fixed point measure of each operator M (l) µ is given by Φ l (µ); that is we have ∀l ≥ 0 ∀µ ∈ P(S (l-1) ) Φ l (µ)(dy) = Φ l (µ)(dx) M (l) µ (x, y)

For l = 0, we use the convention M (0) µ = M (0) and Φ 0 (µ) = π (0) . The initial distributions of the i-MCMC algorithm will be defined in terms of a series of probability measures ν l on S (l) , with l ≥ 0. Now, we are in position to defined inductively the i-MCMC stochastic algorithm. At level l = 0, we let X (0) := (X (0) n ) n≥0 be a Markov chain on S (0) with initial distribution ν 0 and elementary Markov transitions M (0) . For every k ≥ 1, given a realization of the chain X (k-1) := (X (k-1) n ) n≥0 , the k-th level chain X of the chain at level (k -1); that is we have that

P(X (k) n+1 ∈ dx | X (k-1) , X (k) n ) = M (k) η (k-1) n (X k n , dx) (1.2)
with

η (k-1) n := 1 n + 1 n p=0 δ X (k-1) p
The rationale behind this is that the k-th level chain X (k) n behaves asymptotically as a Markov chain with time homogeneous elementary transitions M (k) π (k-1) , as long as η (k-1) n is a good approximation of π (k-1) . For a detailed discussion on the origins and the application model areas of these i-MCMC models we refer the reader to the pair of recent articles [2, 3]. For sufficiently regular models, we proved in [2, 3] that the occupations measures η n (f ) = π (l) , almost surely and for any bounded measurable function f on S (l) . The pair of articles [3, 4] also provide a collection of non asymptotic L rmean error estimates and exponential deviations inequalities. The fluctuation analysis of η (l) n around the limiting measure π (l) has been initiated in [4] in the special case where M (k) µ (x k , .) = Φ k (µ). In this situation, we observe that the k-th level chain X (k) is a collection of independent random variables with distributions ν k for n = 0 and Φ k (η

(k-1)
n-1 ) for n ≥ 1. The fluctuation analysis of this model is essentially based on the study of sequences of conditionally independent random variables X (k) , whose distributions depend in a nonlinear way on the occupation measure of a similar series of random states X (k-1) on the lower level set, and so on for lower indexed levels.

The present article studies the fluctuations of the occupation measures associated with the abstract general class of i-MCMC models discussed above. The theoretical analysis of this type of models is much more involved than the one of traditional Markov chains since it relies on the stability properties of sophisticated Markov chains with elementary transitions that depend in a nonlinear way on the occupation measure of the chain.

The rest of the paper is organized as follows.

The main notation and convention used in this work are presented in a short preliminary section, section 1.2. The main result of the article is presented in full details in section 2. The regularity conditions are summarized in section 2.1. In section 2.2 we state a functional central limit theorem in terms of the semigroup associated with a first order expansion of the one step evolution mappings Φ l of the limiting measure-valued equation (1.1). Section 4 addresses the fluctuation analysis of an abstract class of time inhomogeneous Markov chains. In section 4.2, we present a preliminary resolvent analysis to estimate the regularity properties of resolvent operators and invariant measure type mappings. In section 4.3, we apply these results to study the local fluctuations of a class of weighted occupation measures associated with these models. Section 5 addresses the fluctuation analysis of local interaction random fields associated with the i-MCMC model. The proof of the main theorem presented in section 2.2 is a direct consequence of a fluctuation theorem for local interaction random fields, and it is given at the end of section 5.1.

Notation and conventions

We denote respectively by M(E), M 0 (E), P(E), and B(E), the set of all finite signed measures on some measurable space (E, E), the convex subset of measures with null mass, the subset of all probability measures, and finally the Banach space of all bounded and measurable functions f on E equipped with the uniform norm f = sup x∈E |f (x)|. We also denote by B 1 (E) ⊂ B(E) the unit ball of functions f ∈ B(E) with f ≤ 1, and by Osc 1 (E), the convex set of Emeasurable functions f with oscillations less than one; that is,

osc(f ) = sup {|f (x) -f (y)| ; x, y ∈ E} ≤ 1
We let µ(f ) = µ(dx) f (x), be the Lebesgue integral of a function f ∈ B(E), with respect to a measure µ ∈ M(E). We slight abuse the notation, and sometimes we denote by µ(A) = µ(1 A ) the measure of a measurable subset A ∈ E.

We recall that a bounded integral operator M from a measurable space (E, E) into an auxiliary measurable space (F, F) into itself, is an operator f → M (f ) from B(F ) into B(E), so that the functions

M (f )(x) = F M (x, dy) f (y) ∈ R
are E-measurable and bounded, for any f ∈ B(F ). By Fubini's theorem, we recall that a bounded integral operator M from a measurable space (E, E) into an auxiliary measurable space (F, F) also generates a dual operator µ → µM from M(E) into M(F ) defined by (µM )(f ) := µ(M (f )).

We denote by M := sup f ∈B1(E) M (f ) the norm of the operator f → M (f ) and we equip the Banach space M(E) with the corresponding total variation norm µ = sup f ∈B1(E) |µ(f )|. We let β(M ) be the Dobrushin coefficient of a bounded integral operator M defined by the following formula

β(M ) := sup {osc(M (f )) ; f ∈ Osc 1 (F )}
When M has a constant mass, that is M (1)(x) = M (1)(y) for any (x, y) ∈ E 2 , the operator µ → µM maps M 0 (E) into M 0 (F ), and β(M ) coincides with the norm of this operator. We equip the sets of distribution flows M(E) N with the uniform total variation distance defined by

∀η = (η n ) n≥0 , µ = (µ n ) n≥0 ∈ M(E) N η -µ := sup n≥0 η n -µ n
We extend a given bounded integral operator µ ∈ M(E) → µM ∈ M(F ) into an mapping

η = (η n ) n≥0 ∈ M(E) N → ηD = (η n M ) n≥0 ∈ M(F ) N
Sometimes, we slightly abuse the notation and we denote by ν instead of (ν) n≥0 the constant distribution flows equal to a given measure ν ∈ P(E).

For any R d -valued function f = (f i ) 1≤i≤d ∈ B(F ) d , any integral operator M from E into F , and any µ ∈ M(F ), we will slightly abuse the notation, and we write M (f ) and µ(f ) the R d -valued function and the point in R d given respectively by

M (f ) := M (f 1 ), . . . , M (f d )
and µ(f

) := µ(f 1 ), . . . , µ(f d )
We also simplify the notation and sometimes we write

M [(f 1 -M (f 1 )) (f 2 -M (f 2 ))](x) instead of M [(f 1 -M (f 1 )(x)) (f 2 -M (f 2 )(x))](x) = M (f 1 f 2 )(x) -M (f 1 )(x) M (f 2 )(x)
Unless otherwise is stated, we denote by c(k), k ∈ N, a constant whose value may vary from line to line but only depends on the parameter k. For any pair of integers 0 ≤ m ≤ n, we denote by (n) m := n!/(n -m)! the number of one to one mappings from {1, . . . , m} into {1, . . . , n}. Finally, we shall use the conventions ∅ = 0 and ∅ = 1.

RR n°0123456789

inria-00227536, version 5 -5 Feb 2008

2 Statement of some results

Regularity conditions

Our first regularity condition is a first order weak regularity condition on the mappings Φ l governing the measure-valued equation (1.1). Namely, we further assume that, for any l ≥ 0, the mappings Φ l+1 : P(S (l) ) → P(S (l+1) ) satisfy the following first order local decomposition

Φ l+1 (µ) -Φ l+1 (π (l) ) = (µ -π (l) )D l+1 + Ξ l (µ, π (l) ) (2.1)
In the above display, D l+1 : B(S (l+1) ) → B(S (l) ) stands for a bounded integral operator from S (l) into S (l+1) (that may depend on the measure π (l) ); and Ξ l (µ, η) stands for a remainder signed measure on S (l+1) indexed by the set of probability measures µ, η ∈ P(S (l) ). We further require that

|Ξ l (µ, η)(f )| ≤ (µ -η) ⊗2 (g) Ξ l (f, dg) (2.2)
for some integral operator Ξ l from B(S (l+1) ) into the set T 2 (S (l) ) of all tensor product functions g = i∈I a i (h

1 i ⊗ h 2 i ), with I ⊂ N, (h 1 i , h 2 i ) i∈I ∈ (B(S (l) ) 2 ) I , and a sequence of numbers (a i ) i∈I ∈ R I such that |g| := i∈I |a i | h 1 i h 2 i < ∞ and χ l := sup f ∈B1(S (l) ) |g| Ξ l (f, dg) < ∞
Our second regularity condition concerns the integral operators M β((

M (l) η ) n l ) < 1 (2.3) and M (l) µ -M (l) ν (f ) ≤ |[µ -ν] (g)| Γ l,µ (f, dg) (2.4)
for some collection of bounded integral operators Γ l,µ from B(S (l) ) into B(S (l-1) ), and indexed by the set of measures µ ∈ P(S (l-1) ) with sup

µ∈P(S (l-1) ) Γ l,µ (f, dg) g ≤ Λ l f and Λ l < ∞
We end this section with some comments on this set of conditions. Firstly, we mention that the regularity condition (2.1) is a first order refinement of a Lipschitz type condition we used in [2, 3] to derive a series of L p -mean error bounds and exponential inequalities. This condition has been introduced in [4] for studying the fluctuations of the elementary i-MCMC model when M µ . This conditions is less stringent than the one we used in [3]. In the latter referenced article, we assumed that (2.4) holds for some operators Γ l,µ = Γ l that do not depend on µ.

Therefore, most of the asymptotic results presented in [3] do not apply directly in the present context. Nevertheless, it can be checked that the inductive proof of the L p -mean error bounds presented in theorem 5.2 hold true under the weaker condition (2.4); thus, for every l ≥ 0 and any function f ∈ B(S (l) ), we also know that η

(l) n (f ) converge almost surely to π (l) (f ) as n → ∞.
The main advantage of this set of conditions is that it stable under a state space enlargement, so that the asymptotic analysis of such models, including the functional central limit theorem presented in the next section, applies directly without further work to i-MCMC models on path spaces. We shall return to this key property in section 3.

A functional central limit theorem

To describe precisely the fluctuations of the empirical measures η (l) n around their limiting value π (l) , we need a few notation. We denote by D k,l with 0 ≤ k ≤ l the semigroup associated with the bounded integral operators D k introduced in (2.1). More formally, we have

∀1 ≤ k ≤ l D k,l = D k D k+1 . . . D l
For k > l, we use the convention D k,l = Id, the identity operator. We refer the reader to [4] for an explicit functional representation of these semigroups for a class of Feynman-Kac models arising in physics, biology and information engineering sciences.

In this notation, the functional central limit theorem describing the fluctuations of the i-MCMC model around the solution of the equation (1.1) is stated as follows.

Theorem 2.1 For every k ≥ 0, the sequence of random fields (U

(k) n ) n≥0 on B(S (k) ) defined below U (k) n := √ n η (k) n -π (k)
converges in law, as n tends to infinity and in the sense of finite dimensional distributions, to a sequence of Gaussian random fields U (k) on B(S (k) ) given by the following formula

U (k) := 0≤l≤k (2l)! l! V (k-l) D (k-l)+1,k
In the above display, V (l) l≥0 stands for a collection of independent and centered Gaussian fields with a variance function given by

E V (l) (f ) 2 = π (l) (f -π (l) (f )) 2 + 2 n≥1 π (l) (f -π (l) (f )) M (l) π (l-1) n (f -π (l) (f )) (2.5)
We end this section with some comments on this functional central limit theorem.

RR n°0123456789

inria-00227536, version 5 -5 Feb 2008

Firstly, in the special case where M (l) µ (x l , .) = Φ l (µ) for all l ≥ 1, the variance function (2.5) reduces to the one we obtained in [4]; that is we have

E V (l) (f ) 2 = π (l) (f -π (l) (f )) 2
Loosely speaking, this special class of i-MCMC model behaves as a sequence of independent random variables with distributions Φ l (η

(l-1) n
) given by the local invariant measures of Markov chain Monte Carlo model with general transitions

M (l) η (l-1) n (x l , .).
In these interpretations, the additional terms in the right hand side (2.5) reflects the fluctuations of this MCMC model around their limiting invariant measure. Although the choice M (l) µ (x l , .) = Φ l (µ) seems in some sense optimal, it is difficult to compare the fluctuation variances of these models.

Our second remark is that the functional central limit theorem stated above applies directly without further work to the analysis of the fluctuations of the occupation measures of (X (k) n ) 0≤k≤l around the limiting tensor product measure ⊗ 0≤k≤l π (k) , for any time horizon l ≥ 0. We shall discuss this stability property in the next section devoted to i-MCMC models on path spaces.

Since preparing the original version of this paper, the authors' attention have been drawn to a recent work of Yves Atchadé [1] which analyzes the local fluctuations of some related algorithms; namely the importance-resampling MCMC algorithm and a version of the equi-energy sampler. In [1], the author provides a pair of local central limit theorems for quantities associated to some random measures π (k) n which converge almost surely to π (k) as n tends to infinity, and given by the following formula

1 √ n + 1 n p=0 f k (X (k) p ) -π (k) p (f k )
In our notational system, it is shown in [1] that the above quantities converge in law to a centered Gaussian random variable V (k) (f k ) that does not depend on the first order semigroup D k,l , nor on the fluctuations of the i-MCMC model on lower indexed levels. As it will become transparent in the further developments of the present article, this local fluctuation result is only true if we choose the random sequence π

(k) p := Φ k (η (k-1)
p-1 ), and it does not reflect the complete fluctuations of the i-MCMC model around the limiting measures π (l) .

Path space i-MCMC models

In this section, we fix a final time horizon l of the equation (1.1) and we design an i-MCMC interpretation of the path space model given by ∀n ≥ 0

X [l] n := (X (0) n , . . . , X (l) n ) ∈ S [l] := 0≤k≤l S (k)
For every 0 ≤ k ≤ l we denote by η (k) ∈ P(S (k) ) the image measure of a measure η ∈ P(S [l] ) on the k-th coordinate level set S (k) of the product space k) . In this notation, it is readily checked that X

S [l] := 0≤k≤l S (
[l]

n is an S [l]valued self-interacting Markov chain with elementary transitions defined by

P(X [l] n+1 ∈ dx | (X [l-1] p ) 0≤p≤n , X [l] n ) = M [l] η [l-1] n (X [l] n , dx) (3.1)
with the occupation measures η

[l-1] n
and the collection of transitions

M [l] η [l-1] n
defined by the following synthetic formulae

η [l-1] n := 1 n + 1 n p=0 δ X [l-1] p and M [l] η [l-1] n (X [l] n , dx) = 0≤k≤l M (k) η [l],(k-1) n (X (k) n , dx l )
In the above displayed formulae we have use the convention

M (0) η (-1) n = M (0)
for l = 0, and dx := dx 0 × . . . × dx m stands for an infinitesimal neighborhood of a generic point in the product space x := (x 0 , . . . , x l ) ∈ S [l] . We leave the reader to check that (3.1) coincides with the i-MCMC model associated with the limiting evolution equation

π [l] = Φ [l] (π [l-1] ) with π [l] := π (0) ⊗ . . . ⊗ π (l)
and the invariant measure mapping

Φ [l] : µ ∈ P(S [l-1] ) → Φ [l] (µ) := π (0) ⊗ Φ 1 (µ (0) ) ⊗ . . . ⊗ Φ l (µ (l-1) ) ∈ P(S [l] )
To describe the main result of this section, we need to introduce additional notation. For any 0 ≤ k 1 ≤ k 2 , we set

S [k1,k2] := k1≤k≤k2 S (k) and π [k1,k2] := ⊗ k1≤k≤k2 π (k) ∈ P(S [k1,k2] )
For any 0 ≤ k < l, any pair (µ 1 , µ 2 ) ∈ P(S [0,k] ) × P(S [k+2,l+1] ) and any integral operator D from S (k) into S (k+1) , we denote by (µ

1 ⊗ D ⊗ µ 2 ) from S [l] into S [l+1] (µ 1 ⊗ D ⊗ µ 2 )((x 1 , x 2 , x 3 ), d(y 1 , y 2 , y 3 )) = µ 1 (dy 1 ) D(x 1 , dy 2 ) µ 2 (dy 3 )
where d(y 1 , y 2 , y 3 ) = dy 1 × dy 2 × dy 3 stands for an infinitesimal neighborhood of a generic point (y where 1 ≤ k ≤ l. In addition, the mappings Φ [l+1] satisfy the first order decomposition (2.1) with bounded integral operators D [l+1] from S [l] into S [l+1] given by

1 , y 2 , y 3 ) ∈ S [l+1] = S [0,k] × S (k+1) ⊗ S [k+2,l+1] , and (x 1 , x 2 , x 3 ) a point in the set S [l] = S [0,k-1] × S (k) ⊗ S [k+1,l] .
D [l+1] = π [0,l] ⊗ D l+1 + 0≤k<l π [0,k] ⊗ D k+1 ⊗ π [k+2,l+1] RR n°0123456789

inria-00227536, version 5 -5 Feb 2008

Before presenting the proof of this proposition, we emphasize that this latter directly implies that the functional central limit theorem stated in section 2.2 is also valid for the path space i-MCMC model discussed above. In other words, for every k ≥ 0, the sequence of random fields (U

[k] n ) n≥0 on B(S [k] ) defined by U [k] n := √ n η [k] n -π [k] = 1 √ n + 1 n p=0 δ (X (0 p ,...,X (k) p ) -π (0) ⊗ . . . ⊗ π (k)
converges in law, as n tends to infinity and in the sense of finite dimensional distributions, to a sequence of Gaussian random fields U [k] defined as U (k) by replacing the the semigroups D l1,l2 and the limiting measures π (l) by the corresponding objects on path spaces. Proof of proposition 3.1 : With some obvious notation, we have

∀n ≥ 1 M [l] µ n = 0≤k≤l M (k) µ (k-1)
n

Using the fact that

0≤k≤l µ k - 0≤k≤l η k ≤ 0≤k≤l µ k -η k
for any sequence of probability measures µ k , η k ∈ P(S (k) ), with 0 ≤ k ≤ l, we prove that

β M [l] µ n ≤ 0≤k≤l β M (k) µ (k-1) n ≤ 1 l + 1 0≤k≤l m k (n k ) < 1 as soon as n ≥ n [l] := (∨ 0≤k≤l n k ) × 1 + log (l + 1) ∧ 0≤k≤l log (1/m k (n k ))
We prove the pair of regularity conditions (2.1) and (2.4) by induction on the parameter l. We use the notation

D [l] , Ξ [l] , Ξ [l] m [l] , n [l]
and Γ [l],µ the corresponding objects introduced in the statement of conditions (2.1), (2.3), and (2.4). The results are clearly true for m = 0 with

Φ [0] (µ) := π (0) and M [0] η [-1] n := M (0)
In this case, we readily find that

D [0] = D 0 = 0, Ξ [l] = 0, m [0] = m 0 , n [0] = n 0 and Γ [0],µ = Γ 0,µ = 0
Assume now that the result has been proved at some rank l. For any measure µ on S [l+1] = S [l] × S (l+1) , we denote by µ [l] and µ (l+1) its image measures on S [l] and S (l) . In this notation we have

Φ [l+1] (µ) = Φ [l] (µ [l] ) ⊗ Φ l+1 (µ (l) )
and After some elementary computations, using the decomposition

M [l+1] µ ((u, v), d(x, y)) = M [l] µ [l] (u, dx) × M (l) µ (l) (v,
Φ [l+1] (µ) -Φ [l+1] (π (l) ) = Φ [l] (µ [l-1] ) ⊗ Φ l+1 (µ (l) ) -Φ [l] (π [l-1] ) ⊗ Φ l+1 (π (l) ) = Φ [l] (π [l-1] ) ⊗ Φ l+1 (µ (l) ) -Φ l+1 (π (l) ) + Φ [l] (µ [l-1] ) -Φ [l] (π [l-1] ) ⊗ Φ l+1 (π (l) ) + Φ [l] (µ [l-1] ) -Φ [l] (π [l-1] ) ⊗ Φ l+1 (µ (l) ) -Φ l+1 (π (l) )
we find that the first order condition (2.1) is satisfied with an integral operator l+1] defined for any f ∈ by

D [l+1] from S [l] into S [
D [l+1] (u, d(x, y)) = Φ [l] (π [l-1] )(dx) D (l+1) (u, dy) + D [l] (u, dx) Φ l+1 (π (l) )(dy) = π [l] (dx) D (l+1) (u, dy) + D [l] (u, dx) π (l+1) (dy)
Condition (2.4) is proved using the same type of arguments. This ends the proof of the proposition.

On the fluctuations of time inhomogeneous Markov chains

We discuss in this section the fluctuations of time inhomogeneous Markov chains with elementary Markov transitions that may depend on some predictable flow of distributions on some possibly different state space.

Description of the models

We consider a collection of Markov transitions M η on some measurable space (S, S) indexed by the set of probability measures η ∈ P(S ), on some possibly different measurable space (S , S ). We further assume that there exists an integer n 0 ≥ 0 such that

m(n 0 ) := sup η∈P(S ) β(M n0 η ) < 1 and we set p(n 0 ) := 2n 0 /(1 -m(n 0 )) (4.1)
We also assume that for any pair of measures (η, µ) ∈ P(S ) 2 we have

|[M µ -M η ] (f )| ≤ |[µ -η] (g)| Γ µ (f, dg) (4.2)
for some collection of bounded integral operator Γ µ from B(S) into B(S ), indexed by the set of measures µ ∈ P(S ) with sup µ∈P(S )

Γ µ (f, dg) g ≤ Λ f for some finite constant Λ < ∞ RR n°0123456789

inria-00227536, version 5 -5 Feb 2008

We associate to the collection of transitions M η an S-valued non homogeneous Markov chain X n with a prescribed initial distribution ν ∈ P(S), and some elementary transitions defined by

∀n ≥ 0 P(X n+1 ∈ dx | F n ) = M ηn (X n , dx) (4.3)
In the above displayed formula, F n is a σ-field that contains the σ-field generated by the random states X p from the origin p = 0 up to the current time horizon p = n. In other words, the randomness included in η n may come from a different stochastic process but the distributions η n are measurable with respect to F n . We further assume that the variations of the flow η n are controlled by some sequence of random variables τ (n) in the sense that For n = 0 we use the convention η -1 = 0, the null measure on S .

∀n ≥ 0 η n -η n-1 ≤ τ (n),

Regularity properties of resolvent operators

The main simplification of conditions (4.1) comes from the fact that M η has an unique invariant measure

Φ(η)M η = Φ(η) ∈ P(S)
In addition, the resolvent operators

P η : f ∈ B(S) → P η (f ) := n≥0 M n η -Φ(η) (f ) ∈ B(S) (4.5)
are well defined absolutely convergent series that satisfy the Poisson equation given by (M η -Id)P η = (Φ(η) -Id) Φ(η)P η = 0

In the rest of the article, we shall use the following definition. Resolvent operators are classical tools for the asymptotic analysis of time inhomogeneous Markov chains. In our context the Markov chain interacts with a flow a probability measures. To analyze the situation where this flow converges to some limiting measure, it is convenient to study the regularity properties of the resolvent operators P η as well as the ones of the invariant measure mapping Φ(η) associated with M η . In addition, for any f ∈ B(S) and any (µ, η) ∈ P(S ) we have the following Lipschitz type inequalities

|[Φ(η) -Φ(µ)] (f )| ≤ |[η -µ] (g)| Υ µ (f, dg) (4.6)
and

[P η -P µ ] (f ) ≤ |[η -µ] (g)| Υ µ (f, dg) (4.7)
where (Υ µ , Υ µ ) is a pair of bounded integral operators from B(S) into B(S ) indexed by the set of measures µ ∈ P(S ) such that

g Υ µ (f, dg) ≤ p(n 0 ) Λ f and g Υ µ (f, dg) ≤ p(n 0 )(1 + p(n 0 )) Λ f Proof:
The first assertion is proved in [3]. For completeness, it is only sketched here. We simply use the fact that

P η (f )(x) = n≥0 M n η (f )(x) -M n η (f )(y) Φ(η)(dy) to check that P η (f ) ≤ n≥0 osc(M n η (f ))
and

P η (f ) ≤   n≥0 β(M n η )   osc(f ) ⇒ P η ≤ 2 p≥1 n0-1 r=0 β(M (p-1)n0+r η ) ≤ 2n0 1-β(M n 0 η )
The end of the proof of the first assertion is now immediate. The proof of (4.6) is based on the following decomposition

[Φ(η) -Φ(µ)] (f ) = {[Φ(η) -Φ(µ)] M µ + Φ(η) [M η -M µ ]} (f µ ) with f µ := (f -Φ(µ)(f ))
. Under our regularity conditions on the integral operators M µ , we find that

| [Φ(η) -Φ(µ)] (f )| ≤ | [Φ(η) -Φ(µ)] M µ (f µ )| + [M η -M µ ] (f µ ) ≤ | [Φ(η) -Φ(µ)] M µ (f µ )| + |[µ -η] (g)| Γ µ (f µ , dg) RR n°0123456789
inria-00227536, version 5 -5 Feb 2008 Definition 4.3 We let W be the set of non negative and non increasing weight array functions w = (w n (p)) 0≤p≤n,0≤n , satisfying the following conditions

lim n→∞ ∨ 0≤p≤n w n (p) = 0 with ∀ ∈ [0, 1] ( ) := lim n→∞ 0≤p≤ n w 2 n (p) < ∞
and some scaling function such that lim ( 0, 1 )→(0-,1+) ( ( 0), ( 1)) = (0, 1). We observe that the traditional and constant fluctuation rates sequences w n (p) = 1/ √ n belong to W, with the identity function ( ) = .

Definition 4. 4 We associate to a given weight array function w ∈ W the mapping

W : η ∈ M(S) N → W (η) = (W n (η)) n≥0 ∈ M(S) N
defined for any flow of measures η = (η n ) n≥0 ∈ P(S), and any n ≥ 0, by the weighted measures

W n (η) = 0≤p≤n w n (p) η p
The next proposition presents a pivotal decomposition formula of the weighted occupation measures in terms of a particular martingale on fixed time horizon with a negligible remainder bias term. For n = 0, we use the convention Φ (η -1 ) = ν 0 , so that ζ 0 = [δ X0 -ν 0 ]. For any weight array function w ∈ W, the weighted measures W n (ζ) satisfy the following decomposition

W n (ζ)(f ) = 0≤p≤n w n (p) ∆M p+1 (f ) + L n (f ) (4.8)
with the martingale increments

∆M p+1 (f ) = M p+1 (f ) -M p (f ) := P ηp-1 (f )(X p+1 ) -M ηp P ηp-1 (f )(X p )
and a remainder signed measure L n such that

L n ≤ w n (0) (1 + p(n 0 )τ (n)) (2 + p(n 0 )) Λ Proof:
We let P ηn-1 be the integral operator solution of the Poisson equation associated with the Markov transition M ηn-1 , with a fixed point measure Φ (η n-1 ). By construction, we have

ζ n (f ) = [f (X n ) -Φ (η n-1 ) (f )] = P ηn-1 (f )(X n ) -M ηn-1 (P ηn-1 (f ))(X n )
For n = 0, we use the convention P η-1 = Id and M η-1 = ν 0 . The proof of (4.8) is based on the following decomposition

ζ n (f ) = A n (f ) + B n (f ) + C n (f ) + ∆M n+1 (f ) RR n°0123456789
inria-00227536, version 5 -5 Feb 2008

with the random processes A n (f ), B n (f ) and C n (f ) defined below

A n (f ) := P ηn -P ηn-1 (f )(X n+1 ) B n (f ) := P ηn-1 (f )(X n ) -P ηn (f )(X n+1 ) C n (f ) := M ηn -M ηn-1 P ηn-1 (f )(X n )
Using the Lipschitz inequality (4.7) presented in proposition 4.2, we prove that

|A n (f )| ≤ P ηn -P ηn-1 (f ) ≤ |[η n -η n-1 ] (g)| Υ ηn (f, dg) ≤ τ (n) p(n 0 )(1 + p(n 0 )) Λ f
In addition, using the Lipschitz regularity condition (4.2), we also obtain

|C n (f )| ≤ M ηn -M ηn-1 (P ηn-1 (f )) ≤ |[η n -η n-1 ] (g)| Γ ηn (P ηn-1 (f ), dg) ≤ τ (n) Λ P ηn-1 f ≤ τ (n) Λ p(n 0 ) f
By definition of the weighted measure W n (ζ), we have

W n (ζ) := 0≤p≤n w n (p)ζ p (f ) = 0≤p≤n w n (p) ∆M p+1 (f ) + 0≤p≤n w n (p)(A p (f ) + B p (f ) + C p (f )) (4.9)
From previous calculations, we have

0≤p≤n w n (p)(A p (f ) + C p (f )) ≤ w n (0) τ (n) Λ p(n 0 )(2 + p(n 0 )) f
Finally, we use the following decomposition

0≤p≤n w n (p)B p (f ) = 0≤p≤n w n (p) P ηp-1 (f )(X p ) -w n (p + 1)P ηp (f )(X p+1 ) + 0≤p≤n [w n (p + 1) -w n (p)] P ηp (f )(X p+1 )
with the convention w n (n + 1) = 0. This readily implies that

0≤p≤n w n (p)B p (f ) ≤ 2 w n (0) f + p(n 0 ) f 0≤p≤n [w n (p) -w n (p + 1)] = (2 + p(n 0 )) f w n (0)
The end of the proof is now a direct consequence of formula (4.9). Now, we are in position to state and to prove the main result of this section.

INRIA

Theorem 4. 6 We suppose there exists some measure η, and some m ≥ 1 such that

∀f ∈ B 1 (S ) E(|η n (f ) -η(f )| m ) ≤ m (n) with n≥0 m (n) < ∞
We let V n := W n (ζ) be the sequence of random fields on B(S) associated with a given weight array function w ∈ W and defined in (4.8). We suppose that w ∈ W is chosen so that w n (0)τ (n) tends to 0 as as n → ∞. In this situation,

V n converges in law, as n → ∞ to a Gaussian random field V on B(S) such that ∀(f, g) ∈ B(S) 2 E(V (f )V (g)) = Φ(η) [C η (f, g)] with the local covariance functions C η (f, g) := M η [(P η (f ) -M η P η (f )) (P η (g) -M η P η (g))]
Proof: Using proposition 4.5, it is clearly sufficient to prove that the random fields

W n (ζ) := 0≤p≤n w n (p) ∆M p+1 (4.10)
converge in law to the Gaussian random field V as n → ∞. To use the Lindeberg central limit theorem for triangular arrays of R d -valued random variables, we let f = (f i ) 1≤i≤d ∈ B(S) d be a collection of d-valued functions; and we consider the

R d -valued random variables W n (ζ)(f ) = (W n (ζ)(f i )) 1≤i≤d .
We further denote by F p the σ-field generated by the random variables X q for any q ≤ p. By construction, for any functions f and g ∈ B(S) and for every 0 ≤ p ≤ n we find that

E(w n (p) ∆M p+1 (f ) | F p ) = 0 E(w n (p) 2 ∆M p+1 (f )∆M p+1 (g) | F p ) = w n (p) 2 C p (f, g)(X p )
with the local covariance function

C p (f, g) := M ηp P ηp-1 (f ) -M ηp P ηp-1 (f k ) P ηp-1 (g) -M ηp P ηp-1 (g)
Using proposition 4.2, after some tedious but elementary calculations we find that

C p (f, g) -C η (f, g) ≤ c(η) |[η p-1 -η] (h)| Υ 1 η ((f, g), dh) + |[η p -η] (h)| Υ 2 η ((f, g), dh) with pair of bounded integral operator Υ i η , i = 1, 2, from B(S) 2 into B(S ) such that h Υ i η ((f, g), dh) ≤ c(η) f g
In the above displayed formula, c(η) < ∞ stands for a finite constant whose value only depends on the measure η. Under our assumptions, the following almost sure convergence result readily follows

lim p→∞ C p (f, g) -C µ (f, g) = 0 (4.11) RR n°0123456789
inria-00227536, version 5 -5 Feb 2008

We let π n be the distribution flow of the random variables X n . We proved in [3] that the quantities

1 n + 1 0≤p≤n (h(X p ) -Φ(η p )(h))
converge almost surely to 0 for any function h ∈ B(S) as n → ∞. On the other hand, we also have the almost sure convergence result

lim n→∞ 1 n + 1 0≤p≤n Φ(η p )(h) = Φ(η)(h)
This also implies that

lim n→∞ E   1 n + 1 0≤p≤n h(X p )   = lim n→∞ 1 n + 1 0≤p≤n π p (h) = Φ(η)(h)
from which we conclude that lim n→∞ π n (h) = Φ(η)(h). On the other hand, using the Burkholder-Davis-Gundy type inequalities for martingales, we find that for any weight array function (v n (p)) 0≤p≤n and any m ≥ 1, we have

E   0≤p≤n v n (p)(h(X p ) -π p (h)) m   1 m ≤ a(m)   0≤p≤n v n (p) 2   1 2 osc(h)
for some constants a(m) whose values only depend on the parameter m. Thus, if we set v n (p) = w n (p) [START_REF] Brockwell | Sequentially interacting Markov chain Monte Carlo[END_REF] , we find that

E 0≤p≤n w 2 n (p)(h(X p ) -π p (h)) m 1 m ≤ a(m) 0≤p≤n w n (p) 4 1 2 osc(h) ≤ a(m) w n (0) 0≤p≤n w n (p) 2 1 2 osc(h)
Under our assumptions on the weight functions w, if we take h = C η (f, g) then by (4.11) we obtain the following almost sure convergence result

lim n→∞ n p=0 w n (p) 2 C p (f, g)(X p ) = lim n→∞ n p=0 w n (p) 2 C η (f, g)(X p ) = lim n→∞ n p=0 w n (p) 2 π p (C η (f, g)) = Φ(η)(C η (f, g))
Therefore, we also have the almost sure convergence result

lim n→∞ n p=0 w n (p) 2 E(∆M p+1 (f )∆M p+1 (g) | F p ) = Φ(η)(C η (f, g))
Since we have ∨ 0≤p≤n w n (p) = w n (0) → 0, as n → ∞, the Lindeberg condition is satisfied and we conclude that the sequence of random fields W n (ζ) defined in (4.10) converges in law to the Gaussian random field V as n → ∞. This ends the proof of the theorem.

We end this section with an alternative and simpler representation of the covariance function of the random field V presented in proposition 4.6. We have

C η (f, f )(x) = M η (x, dy) [P η (f )(y) -M η (P η (f ))(x)] 2
Using the decomposition

P η (f )(y) -M η (P η (f ))(x) = [P η (f )(y) -P η (f )(x)] + [P η (f )(x) -M η (P η (f ))(x)] = [P η (f )(y) -P η (f )(x)] + [f (x) -Φ(η)(f )]
and the fact that

M η (x, dy) [P η (f )(y) -P η (f )(x)] = [M η (P η (f ))(x) -P η (f )(x)] = -[f (x) -Φ(η)(f )]
we readily prove the formula

C η (f, f )(x) = M η (x, dy) [P η (f )(y) -P η (f )(x)] 2 -[f (x) -Φ(η)(f )] 2
On the other hand, recalling that Φ(η) = Φ(η)M η and using again the Poisson equation we also have

Φ(η)(dx)M η (x, dy) [P η (f )(y) -P η (f )(x)] 2 = 2 Φ(η) [P η (f ) (f -Φ(η)(f ))] and 2 Φ(η) [P η (f ) (f -Φ(η)(f ))] = 2 Φ(η) (f -Φ(η)(f )) 2 + 2 n≥1 Φ(η) M n η (f -Φ(η)(f )) (f -Φ(η)(f ))
In summary, we have proved the following proposition.

Proposition 4.7 The limiting covariance function presented in theorem 4.6 is alternatively defined for any function f ∈ B(S) by the following formula

Φ(η) [C η (f, f )] = Φ(η) (f -Φ(η)(f )) 2 + 2 n≥1 Φ(η) (f -Φ(η)(f )) M n η (f -Φ(η)(f ))
∀l ≥ 0 ∀n ≥ 0 δ (l) n := δ X (l) n -Φ l η (l-1) n-1
For n = 0, we use the convention Φ l (η -1 ) = ν l so that δ

(l) 0 = δ X (l) 0 -ν l .
We associate with a sequence of weight array functions (w (l) ) l≥0 ∈ W N , the flow of random fields (W

n (δ (l) )) l≥0 on the sets of functions (B(S (l) )) l≥0 given by

∀l ≥ 0 ∀n ≥ 0 W (l) n (δ (l) ) := 0≤p≤n w (l) n (p) δ (l) p
We observe that the regularity conditions (2.3) and (2.4) ensure that the collection of Markov operators M (l) η and their fixed point measures mappings Φ l (η) satisfy the regularity conditions (4.1) and (4.2) introduced in section 4.1. Also observe that the i-MCMC chain X (l+1) is a time inhomogeneous model of the form (4.3) with a collection of elementary transitions M n associated with the i-MCMC chain at level l. Note that in this situation, condition (4.4) is satisfied with

∀n ≥ 0 η (l) n -η (l) n-1 ≤ τ (l) (n) := 2 n + 1 and we have τ (l) (n) := 0≤p≤n τ (l) (p) = 2 0≤p≤n 1 p + 1 ≤ 2(1 + log (n + 1))
Finally, we recall that for any m ≥ 1 we have

∀l ≥ 0 ∀f ∈ B 1 (S (l) ) E(|η (l) n (f ) -π (l) (f )| m ) 1 m ≤ b(m) c(l) 1 √ n + 1 (5.1)
for a collection of finite constants b(m) whose values only depend on the parameter m (see for instance [3]). Using theorem 4.6, we can prove that the random fields

V (l) n := W (l) n (δ (l) ) (5.2)
associated with a given weight array function w (l) ∈ W where

lim n→∞ (w (l) n (0) log (n)) = 0 converges in law to a Gaussian random field V (l) as n → ∞ such that ∀(f, g) ∈ B(S (l) ) 2 E(V (l) (f )V (l) (g)) = π (l) C (l) (f, g) (5.3)
In the above display, the covariance functions C (l) (f, g) are defined in terms of the resolvent operator P and the fixed point measure Φ l (π (l-1) ) = π (l) with the following formula

C (l) (f, g) := M (l) π (l-1) P (l) π (l-1) (f ) -M (l) π (l-1) P (l) π (l-1) (g) P (l) π (l-1) (g) -M (l) π (l-1) P (l) π (l-1) (g)
The main objective of this section is to prove the following theorem.

INRIA

From the previous discussion, we find that

W (n) l (f i ) = M (n) l (f i ) + L (n) l (f i ) with the F (n) -martingale M (n) l (f i ) given below M (n) l (f i ) := 0≤k≤l ∆M (n) k (f i ) with ∆M (n) k (f i ) := 0≤p≤n w (k) n (p) ∆M (k) p+1 (f )
(5.5) and the remainder bias type measure L

(n) l = 0≤k≤l L (k) n , such that lim n→∞ L (n) l = 0
The theorem 5.2 is now a direct consequence of the following proposition (see for instance the arguments used in section 4.2 in [4]).

Proposition 5.3 The sequence of martingales M (n) l (f ) defined in (5.5) converges in law, as n → ∞, to an R d -valued and Gaussian martingale M l (f ) = (M l (f i )) 1≤i≤d such that for any l ≥ 0 and any pair of indexes

1 ≤ i, j ≤ d M(f i ), M(f j ) l = 0≤k≤l π (k) C (k) (f i , f j )
with the local covariance functions π (k) C (k) (f i , f j ) defined in (5.3).

Proof:

The proof of the proposition follows the same lines of arguments as those we used in the proof of theorem 4.6. Firstly, we consider the decomposition

M (n) l (f i ) = l(n+1)+n i=0 V (n) i (f )
where for every 0 ≤ i ≤ l(n + 1) + n, with i = k(n + 1) + p for some 0 ≤ k ≤ l, and 0

≤ p ≤ n V (n) i (f ) := w (k) n (p) ∆M (k) p+1 (f k )
We further denote by G

(n) i the σ-field generated by the pair of random variables (X (k) p , X (k) p+1 ) for any pair of parameters (k, p) such that k(n+1)+p ≤ i. By construction, for any flow of functions f = (f l ) l≥0 and g = (g l ) l≥0 ∈ l≥0 B(S (l) ) and for every 0 ≤ i ≤ l(n + 1) + n, with i = k(n + 1) + p for some 0 ≤ k ≤ l, and 0 ≤ p ≤ n, we find that Using the same arguments as the ones we used in the proof of theorem 4.4 in [4], we conclude that the R d -valued martingale M

E(V (n) i (f ) | G (n) i-1 ) = 0 E(V (n) i (f )V (n) i (g) | G (n) i-1 ) = w (k) n (p) 2 C (k) p (f,
(n) l (f ) converges in law, as n tends to infinity, to a martingale M l (f ) with a predictable bracket given for any air of indexes 1 ≤ j, j ≤ d by

M(f j ), M(f j ) l = C (k) l (f j , f j )
This ends the proof of the proposition.

  chain with initial distribution ν k and with random Markov transitions M (k) η (k-1) n depending on the current occupation measures η (k-1) n

n

  converge to the solution π (l) of the equation (1.1), in the sense that RR n°0123456789 inria-00227536, version 5 -5 Feb 2008 lim n→∞ η (l)

µ

  satisfy the following pair of regularity conditions m l (n l ) := sup η∈P(S (l-1) )

µ

  (x, .) = Φ l (µ). The second regularity condition (2.4) is a local Lipschitz type continuity condition on the collection of operators M (l)

Proposition 3 . 1

 31 For any l ≥ 0, the mappings Φ [l] and the collection of Markov transitions M [l] µ satisfy the set of regularity conditions (2.1), (2.3), and (2.4) as long as the corresponding conditions are met for the marginal mappings Φ k and the transitions M (k) µ

  dy) INRIA inria-00227536, version 5 -5 Feb 2008 for any (u, v) ∈ S [l+1] = (S [l] × S (l) ), where d(x, y) = dx × dy stand for an infinitesimal neighborhood of a generic point (x, y) ∈ S [l+1] = (S [l] × S (l) ).

  and we set τ (n) := 0≤p≤n τ (p) (4.4)

Definition 4 . 1

 41 The integral operator P η defined by the series (4.5) is called the resolvent operator associated to the Markov transition M η and the fixed point measure Φ(η).

Proposition 4 . 2

 42 Under the regularity conditions (4.1) and (4.2), we have sup η∈P(S ) P η ≤ p(n 0 )

Proposition 4 . 5

 45 We consider the flow of random measures ζ := (ζ n ) n∈N ∈ M(S) N defined for any n ≥ 0 by the following formula ∀n ≥ 0 ζ n = [δ Xn -Φ (η n-1 )]

5 A fluctuation theorem for local interaction fields 5 . 1 IntroductionDefinition 5 . 1

 5151 This section presents the fluctuation analysis of a class of weighted random fields associated with the i-MCMC model. Following the local fluctuation analysis for time inhomogenous Markov chain models presented in section 4.3, we introduce the following weighted random fields models.RR n°0123456789inria-00227536, version5 -5 Feb 2008 We consider the flow of random measures

  the flow of occupation measures η (l)

1 (n + 1 0≤p≤nπ 1 2 1 2 2 C

 11112 g k ) RR n°0123456789inria-00227536, version5 -5 Feb 2008 Under our Lipschitz regularity conditions (2.3) and (2.4), the proposition 4.2 applies to the fixed point mappings Φ k and the resolvent operators P (k) µ . As in the proof of theorem 4.6, after some tedious but elementary calculations we obtainC (k) p (f, g) -C (k) (f, g) ≤ c(k) η (k-1) p-1 -π (k-1) (h) Υ (k),1 π k ,π (k-1) ((f k , g k ), dh) + η (k-1) p -π (k-1) (h) Υ (k),2 π k ,π (k-1) ((f k , g k ), dh)whereΥ (k),i π k ,π (k-1) , i = 1, 2, is a pair of bounded integral operator from B(S (k) ) 2 into B(S (k-1) ) such that h Υ (k),i π k ,π (k-1) ((f k , g k ), dh) ≤ c(k) f k g kCombining the generalized Minkowski integral inequality with (5.1) we prove the following almost sure converge resultlim p→∞ C (k) p (f, g) -C (k) (f, g) = 0We let π(k)n be the distribution flow of the random variables X (k)n . One the one hand, since the quantities η (k) n (h) converge almost surely to π (k) (h), as n → ∞, for every k ≥ 0 and any function h ∈ B(S (k) ), we also have thatlim n→∞ E(η (k) n (h)) = lim n→∞ 1 (k) p (h) = π (k) (h) =⇒ lim n→∞ π (k) n (h) = π (k) (h)On the other hand, using the Burkholder-Davis-Gundy type inequalities for martingales, we find that for any weight array function (v n (p)) 0≤p≤n and any m ≥ 1, we have a(m) whose values only depend on the parameter m. Thus, if we set v n (p) = w n (p) 2 , we find that E osc(h) ≤ a(m) w n (0) 0≤p≤n w n (p) 2 osc(h) INRIA inria-00227536, version 5 -5 Feb 2008Under our assumptions on the weight functions w (k) , if we take h = C (k) (f, g) then we obtain the following almost sure convergence resultslim (k) (f, g)(X (k) n (p) 2 π (k) p (C (k) (f, g)) = π (k) (C (k) (f, g))This yields the almost sure convergencelim n→∞ M (n) (f ), M (n) (g) l = C (k) l (f, g) := 0≤k≤l π (k) (C (k) (f, g))

http://www.inria.fr ISSN 0249-6399

This recursion readily implies (4.6) with the integral operator given by Υ µ (f, dg) := n≥0 Γ µ (M n µ (f µ ), dg)

Finally we observe that

Arguing as above, we conclude that

This ends the proof of (4.6). The proof of (4.7) follows the same type of arguments. Firstly, we observe that

To check this formula, we first use the fact that M µ P µ = P µ M µ to prove that

Using the Poisson equation and the fact that P µ (1) = 0 we also have the decomposition P µ (M η -Id)P η = P µ (Φ(η) -Id) = -P µ

Combining these two formulae, we conclude that

This ends the proof of the decomposition given above. It is now easily checked that

The end of the proof follows the same type of arguments as before. This ends the proof of the proposition.

Local fluctuations of weighted occupation measures

This section is concerned with the fluctuation analysis of the occupation measures of the time inhomogeneous Markov chain introduced in (4.3). In section 5, we shall use these results to analyze the fluctuations of i-MCMC models.

The fluctuation analysis of this type of models is related to the fluctuations of weighted occupation measures with respect to some weight array type functions (cf. for instance [4]).

INRIA

Theorem 5.2 We consider a collection of weight array functions (w (l) ) l≥0 ∈ W N where lim n→∞ (w

n (0) log (n)) = 0 for any l ≥ 0. In this situation, the corresponding flow of weighted random fields (V (l) n ) l≥0 defined in (5.2), converges in law, as n tends to infinity and in the sense of finite dimensional distributions, to a sequence of independent and centered Gaussian fields V (l) l≥0 with covariance functions defined in (5.3).

Using this result, the proof of the functional central limit theorem 2.1 follows exactly the same arguments as the one we used in the proof of theorem 2.1 presented in section 6 in the recent article [4].

Proof of theorem 2.1 :

We observe that the time averaged semigroup S k can be rewritten in terms of the following weighted summations s (k) n (q) and s [START_REF] Atchadé | A cautionary tale on the efficiency of some adaptive Monte Carlo schemes[END_REF] n (p) := 1

We also known from proposition 6.1 in [4] that

and for any k ≥ 1, the weight array functions w (k) defined by

belong the set W introduced in definition 4.3, with lim n→∞ (w

On the other hand, using the proposition 5.2 in [4], we have the following multilevel expansion

where

n ) n≥0 is a flow of signed and random measures such that

inria-00227536, version 5 -5 Feb 2008

In the above display, b(m) stand for some constant whose value only depend on the parameter m. This multilevel expansion implies that

with the weighted distribution flow mappings W (k-l) associated with the weight functions w (l+1) and a remainder signed measure Ξ

The proof of theorem 2.1 is now a direct consequence of theorem 5.2. This ends the proof of the theorem.

A martingale limit theorem

This section is mainly concerned with the proof of theorem 5.2. We following the same lines of arguments as the ones we used in section 4.3 devoted to the fluctuations of weighted occupation measures associated with time inhomogeous Markov chain models. Firstly, we fix some notation. For any k ≥ 0 and any µ ∈ P(S (k-1) ), we let P (k) µ be the resolvent operator associated with the Markov transition M (k) µ and the fixed point measure Φ k (µ) ∈ P(S (k) ). We also set

with the pair of parameters (n k , m k ) defined in (2.3).

Using proposition 4.5, we find that the weighted measures W

n (δ (k) ) satisfy the following decomposition

for any f ∈ B(S (k) ) with the martingale increments

and the remainder signed measure L (k)

n which are such that

We consider a sequence of functions f = (f i ) 1≤i≤d , with d ≥ 1, and f i = (f i k ) k≥0 ∈ k≥0 B(S (k) ), and we let W (n) (f ) = (W (n) (f i )) 1≤i≤d be the R dvalued and F (n) -adapted process defined for any l ≥ 0 and any 1 ≤ i ≤ d by W

(n)


Francisco Navas

Vicente Milanés

Mixing V2V-and non-V2V-equipped vehicles in car following

Keywords: Cooperative adaptive cruise control, V2Vequipped vehicles, non-V2V-equipped vehicles, mixed traffic situations, intelligent transportation systems

Cooperative Adaptive Cruise Control (CACC) provides significant traffic flow improvements when a vehicle-tovehicle (V2V) communication link exists with the preceding vehicle, but it degrades to an Adaptive Cruise Control (ACC) when this communication link is no longer available. This degradation occurs even if the information from another V2V-equipped vehicle ahead (different to the preceding one) is still available. This paper presents a novel car-following control system-Advanced Cooperative Adaptive Cruise Control (ACACC)-that benefits from the existing communication with this vehicle ahead in the string, reducing inter-vehicle gap whereas keeping string stability. The proposed control structure provides a hybrid behaviour between two CACC controllers with different time gaps according to the string position of the vehicle with the V2V communication link available. An stable hybrid behavior between both controllers is ensured through the Youla-Kucera parameterization. Simulation and real experiments show the proper behaviour of the designed control algorithm and a good performance compared to existing ACC/CACC controllers.

I. INTRODUCTION

Recently, the International Council on Clean Transportation has estimated that over the next two decades vehicle ownership is expected to increase 7 million only in the European Union (EU) [START_REF]European vehicle market statistics[END_REF]. Therefore, road transport will have to deal with these figures. Among the associated problems, drivers will find more congested roads, resulting in an enormous waste of fuel and productivity together with health problems. According to the European Commission, congestion costs are equivalent to 1% of Gross domestic product (GDP)-in other words more than the EU budget [START_REF] Commision | Road transport. a change of gear[END_REF].

Since building new infrastructure is no longer an appropriate solution, more intelligent and efficient options result of Intelligent Transportation Systems (ITS). Specifically, related to traffic congestion, intelligent longitudinal speed control is a suitable system to improve congestion in highways, through homogeneous speed on the part of the driver and shorter intervehicle distances. Adaptive cruise controller (ACC) [START_REF] Marsden | Towards an understanding of adaptive cruise control[END_REF] is a commercial system already implemented in production vehicles. An ACC system can track the preceding vehicle measuring the actual distance and the ego-vehicle velocity. These inputs allow the system to maintain a selected time gap, calculating the required acceleration or deceleration to reach the desired velocity or to prevent a collision. Recently, research focuses on the cooperative version of the system, so-called cooperative ACC (CACC) [START_REF] Ploeg | Design and experimental evaluation of cooperative adaptive cruise control[END_REF]. Vehicle-to-vehicle F. Navas is with the Research Department, AKKA Technologies, France. The work was done with INRIA Paris, France. {francisco-martin.navas-matos@akka.eu}; V. Milanés is with the Research Department, Renault SAS, 1 Avenue de Golf, 78280

Guyancourt, France {vicente.milanes@renault.com}. (V2V) communication is added to the existing ACC system, improving the traffic flow through the formation of a tighter string of vehicles.

Traffic flow benefits in function of the ACC-equipped vehicles market penetration have been widely studied. Kesting et Al. [START_REF] Kesting | Adaptive cruise control design for active congestion avoidance[END_REF] [START_REF] Kesting | Extending adaptive cruise control to adaptive driving strategies[END_REF] stated that a 25% of market penetration could remove congestion. Recent studies [START_REF] Shladover | Impacts of cooperative adaptive cruise control on freeway traffic flow[END_REF] proposes a different ACC model that provides no traffic flow benefits even with a market penetration of 100%. This model is based on the field test of ACC driven by 16 drivers from the general public in [START_REF] Nowakowski | Cooperative adaptive cruise control: Driver acceptance of following gap settings less than one second[END_REF]. These drivers were encouraged to select the time gap setting that they preferred, resulting in a time-gap close to manual driving. As quantitative example, the highway capacity with every vehicle driven manually is 2050 veh/h, while with ACCequipped vehicles the capacity increases up to 2200 veh/h [START_REF] Werf | Effects of adaptive cruise control systems on highway traffic flow capacity[END_REF].

In contrast, the CACC system on the traffic flow characteristics presents much more optimistic results. Arem et al. [START_REF] Arem | The impact of cooperative adaptive cruise control on traffic-flow characteristics[END_REF] conclude that CACC improves the highway capacity when the penetration rate is greater than 60%, obtaining better results on high traffic volume because of more vehicles participate in the string. The CACC model in [START_REF] Shladover | Impacts of cooperative adaptive cruise control on freeway traffic flow[END_REF] validates these results, showing a maximum lane capacity of 4000 veh/h under 100% CACC-equipped vehicles condition. This model is again based on the chosen time-gap by the general public in the field test in [START_REF] Nowakowski | Cooperative adaptive cruise control: Driver acceptance of following gap settings less than one second[END_REF]. CACC can double the highway capacity in the ideal situation.

It is clear that the market penetration of CACC systems would progressively occur, making necessary to work with mixing traffic situations. As the highway capacity is sensitive to the percentage of CACC-equipped vehicles, it is important to preserve the CACC behaviour when there is no communication with the preceding vehicle, but with another further on the string. For instance, when a vehicle loses communication capability into a string previously formed by CACC-equipped vehicles. As the communication with the preceding vehicle is no longer available, existing literature on the field (see [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] and [START_REF] Ploeg | Graceful degradation of cacc performance subject to unreliable wireless communication[END_REF] for details), degrades the system to a conventional ACC, removing the strong impact of CACC systems on highway capacity. Middle plot of Fig. 1 depicts the situation. As quantitative example, [START_REF] Ploeg | Design and experimental evaluation of cooperative adaptive cruise control[END_REF] showed how the minimum time gap string stable increases from 0.7s to 3.16s when communication is not available.

Recent literature deals with V2V communication problem with the preceding vehicle. Van et Al. [START_REF] Van Nunen | Robust model predictive cooperative adaptive cruise control subject to v2v impairments[END_REF] used a buffer in order to store the communication with the preceding vehicle for a period of time. This buffer is combined with the vehicle model to provide robustness when intermittent communication failures occur. Savaia et Al. [START_REF] Savaia | A receding horizon switching control resilient to communication failures for connected vehicles[END_REF] used a Kalman filter together with the history of preceding vehicle acceleration when intermittent communication failures exist with the preceding vehicle. However, these are useful for short periods of time where communication disappears with the preceding vehicle, otherwise a degradation to ACC needs to be performed for ensuring string stability. CACC controller design when there is no communication with the preceding vehicle, but with another further on, is considered in [START_REF] Gong | Cooperative adaptive cruise control for a platoon of connected and autonomous vehicles considering dynamic information flow topology[END_REF]. A string of three vehicles is employed, giving guidelines in the design of a proportional-derivative (PD) controller and time gap value in order to reduce the effects of communication failures in the string stability of a CACC system. A normal CACC is considered if communication with preceding vehicle is present, otherwise the controller is switched to another one connected to the leader vehicle, avoiding system degradation to an ACC. Nevertheless, stability is not ensured when switching between both controllers, and all vehicles in the string have same dynamics and behavior.

To the best of author knowledge, former CACC working with V2V degradation study either the minimum time gap to guarantee stability according to the V2V lost, switch between CACC and ACC when the V2V link is lost [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] [START_REF] Ploeg | Graceful degradation of cacc performance subject to unreliable wireless communication[END_REF], use prediction for short period of V2V degradation [START_REF] Van Nunen | Robust model predictive cooperative adaptive cruise control subject to v2v impairments[END_REF] [START_REF] Savaia | A receding horizon switching control resilient to communication failures for connected vehicles[END_REF], or give control guidelines (and use communication with a vehicle further on the string) to be able to perform in a mode close to CACC for a longer period of time of degradation, but with all vehicles behaving as CACC [START_REF] Gong | Cooperative adaptive cruise control for a platoon of connected and autonomous vehicles considering dynamic information flow topology[END_REF]. This paper proposes a novel solution when the V2V link with the preceding vehicles is lost without degrading vehicle performance to an ACC. Different preceding vehicle behaviors are considered, it could be either as a CACC-equipped, ACC-equipped or a conventional vehicle. The principle relies on the effect that V2V info from any vehicle ahead is still available so the vehicle can be able to perform in between CACC and ACC responses, leading to shorter inter-vehicular gaps whereas keeping string stability. Besides, Youla-Kucera (YK) parameterization is used in order to ensure stability when adapting the control response between both controllers.

The control algorithm proposed in this work-Advanced CACC (ACACC)-benefits from the V2V information coming from another vehicle ahead in the string instead of degrading to a conventional ACC when there is no communication with the preceding vehicle. Preceding vehicle behavior could be either as a CACC-equipped, ACC-equipped or a conventional vehicle. Since the time gap policy has a significant influence on the highway capacity, ACACC is specifically composed by two CACC controllers with different time gaps. An hybrid behavior between both could not be stable. In [START_REF] Niemann | Switching between multivariable controllers[END_REF], a theoretical example shows how the linear combination of two stabilizing controllers does not result in a stabilizing controller of the system. Different controller interpolation methods have been proposed in the literature to guarantee stability. [START_REF] Shahruz | Design of controllers for linear parametervarying systems by the gain scheduling technique[END_REF] proposed an algorithm for linear interpolation based on eigenplacement of state feedback gains. In [START_REF] Stilwell | Interpolation of observer state feedback controllers for gain scheduling[END_REF] and [START_REF]Stability preserving interpolation methods for the synthesis of gain scheduled controllers[END_REF], approaches for interpolation based on coprime factors and state-space representation are presented. [START_REF] Niemann | An architecture for implementation of multivariable controllers[END_REF] used the YK parameterization of all the controllers to provide controllers' interpolation with guaranteed stability. Later, [START_REF] Hespanha | Switching between stabilizing controllers[END_REF] and [START_REF] Blanchini | A separation principle for linear switching systems and parametrization of all stabilizing controllers[END_REF] guaranteed stability under arbitrary change of the signal responsible of the interpolation over time. Several practical applications of YK parameterization can be found in the project called Plug&Play [START_REF] Stoustrup | Plug & play control: Control technology towards new challenges[END_REF]. YK is chosen among the others to obtain a stable interpolation between the two CACC controllers that compose the proposed ACACC system.

Traffic flow improvement comes from a hybrid behaviour between both while using the existing communication link with the closer ahead communicated vehicle. As the system is not degraded to ACC, the highway capacity previously accomplished by the CACC string of vehicles is better preserved, as well as the shockwave effect is better mitigated. The bottom part of Fig. 1 shows an ACACC scenario; where the red vehicle could be either a conventional driver or an ACCequipped vehicle, both without communication capabilities.

The proposed system has been simulated for highway velocities by using the vehicle model in [START_REF] Milanés | Cooperative adaptive cruise control in real traffic situations[END_REF]. ACACC has been tested under different situations, being the preceding vehicle either a conventional or an ACC-equipped one. For all cases, the system exhibits a good performance, achieving a tighter string of vehicles. These simulations results are later validated on the INRIA experimental platform at low-speeds. ACACC improves traffic flow capacity.

In brief, the main contributions in this paper are the following:

• ACACC benefits for any V2V-equipped vehicle ahead when communication with the preceding one is lost, avoiding a full ACC degradation. Communication unavailability could be either intermittent or permanent. • ACACC is composed by two CACC controllers with different time gaps. YK parameterization is used with these two controllers in order to guarantee a stable interpolation between both. • Scalar factor γ, in charge of interpolation between two controllers in a YK approach, is designed in order to improve traffic flow (reduce inter-vehicle distance, while being string stable). • The approach is adaptable for any type/number of preceding vehicles. Simulation results are provided, being the preceding vehicle either a conventional vehicle emulated by an Intelligent Driver Model (IDM), an ACC-equipped one, or even several ACC-equipped vehicles; all of them without communication capabilities. Simulations results are validated at INRIAs low-speed experimental platform. The rest of the paper is structured as follows. Section II introduces the problem formulation, including the basis of a CACC control algorithm and some initial assumptions. The theory for ensuring a stable interpolation between two general controllers is described in Section III. How this theory is applied to ACACC scenario and how the percentage of each CACC controller is modified through γ for improving highway capacity is presented in Section IV. Section V shows simulation results. The control algorithm is validated on the low-velocity experimental platform in Section VI. Finally, some concluding remarks are given.

II. PROBLEM FORMULATION

A string of n vehicles driving in the same lane, one after the other, is considered. i determines the order of a vehicle inside of the string i ∈ [1, n]. Vehicle 1 is the leader and first vehicle in the string, vehicle i denotes ego-vehicle and vehicle i -1 is the vehicle preceding the ego-vehicle. Leader and preceding vehicle are the same vehicle when n = 2. If n ≥ 3, the number of vehicle between leader 1 and ego vehicle i is n-2. Problem formulation focuses on ego-vehicle situation according to V2V communication availability with vehicle i -1 and the rest of vehicles in the string. V2V communication between vehicle i and vehicle i -1 is defined as C i-1 . Communication between vehicle i and some other vehicle x different from the preceding one is defined as C x with 1 ≤ x < i -1.

In the literature, communication availability with vehicle i -1 determines the use of ACC or CACC controllers in the vehicle i, and therefore its time gap policy h i :

• When C i-1 exists, a regular CACC controller can be used.

Faster responses can be achieved, allowing a tighter string of vehicles by employing a time gap h CACC -i.e. ∃C i-1 → CACC where h i = h CACC . • If C i-1 is no longer available: the system degrades to a conventional ACC algorithm [START_REF] Navas | Using plug&play control for stable acc-cacc system transitions[END_REF]. Consequetly, benefits of CACC systems on highway capacity are removed because a longer time gap h ACC > h CACC needs to be set to ensure string stability-i.e. C i-1 → ACC where h i = h ACC > h CACC . These situations have been extensively studied (see [START_REF] Marsden | Towards an understanding of adaptive cruise control[END_REF] for ACC and [START_REF] Ploeg | Design and experimental evaluation of cooperative adaptive cruise control[END_REF] for CACC for details), but potential benefits of using C x when C i-1 is not available have not been further investigated. Speed oscillations on the string (i.e. the difference between speeds V x and v i-1 ) is limited to 5m/s according to [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF]. Experimental results demonstrated that over this value, drivers disengage ACC system due to its degraded performance.

The objective of the ACACC controller is to enhance traffic flow when there is no communication with vehicle i -1 (preceding vehicle), by taking the information from the closer V2V-equipped vehicle ahead, in that case, vehicle x (that could be either the leader vehicle or another one between preceding and leader). Since V2V communication is always available (no matter from which vehicle comes from), a CACC control structure is proposed (see Sec. II.A for details). ACACC is composed by two CACC controllers with different time gaps: A CACC controller with a short time gap h i = h CACC so-called CACC ST G when C i-1 is available; and a CACC controller with a longer time gap h i = h ACC so-called CACC LT G when C i-1 is lost but C x link is still available. The proposed ACACC algorithm benefits from the YK parameterization to provide a hybrid response between both controllers, assuring control system stability. The regulation between the effect of each controller in the control command is based on a correlation between the speed received via V2V communication and the one detected by the on-board ego-vehicle sensor (i.e. radar). The correct tuning of γ with a maximum speed oscillation of 5m/s, assures string stability and improves traffic flow.

In brief, ACACC works when C i-1 is lost but C x is available, providing a control structure of the form:

i f C i-1 &∃C x ACACC = CACC ST G (1 -γ) +CACC LT G γ; γ ∈ [0, 1] (1)
where

h i ∈ [h CACC , h ACC ].
Notice how the time gap would be lower than degrading the system to a regular ACC.

A. CACC

Figure 2 shows the classical structure for CACC systems (see [START_REF] Ploeg | Design and experimental evaluation of cooperative adaptive cruise control[END_REF] for more details), where G represents the vehicle model, K i a stabilizing controller and V c i-1 the control velocity of vehicle i -1 communicated through C i-1 by means of a feedforward filter F i . The controller should be able to regulate the error e i between relative distance d r,i (X i-1 -X i , as difference between absolute positions) and the reference distance coming from the spacing policy (V i h i , where V i is the ego-vehicle velocity). In previous study [START_REF] Naus | String-stable cacc design and experimental validation: A frequency-domain approach[END_REF], PD controller has provided good performance. PD controller expression is in equation [START_REF] Commision | Road transport. a change of gear[END_REF]. Controller parameters k p,i and k d,i will change for CACC ST G and CACC LT G in the proposed ACACC approach,

K i = k p,i + k d,i s (2)
Highway capacity improves when tighter gap policies can be adopted whereas keeping string stability. String stability is defined as the attenuation of disturbances along the string of vehicles. A sufficient condition for string stability is given in [START_REF] Rajamani | Vehicle dynamics and control[END_REF], which means that the absolute position of each vehicle must not increase as it propagates through the string. Condition is equivalent to equation (3).

X i /X i-1 ∞ ≤ 1 f or i > 1 (3)
According to [START_REF] Naus | String-stable cacc design and experimental validation: A frequency-domain approach[END_REF], for a PD-based CACC system the string stability function results:

X i X i-1 = D i + (1 + h i s)GK i (1 + h i s)(1 + (1 + h i s)GK i ) , f or i > 1 (4)
Under the ideal situation where the communication delay is null (D i = 1), string stability yields:

X i X i-1 = 1 (1 + h i s) , f or i > 1 (5)
Therefore, string stability is guaranteed for any h i > 0. The existence of C i-1 makes vehicle response faster, allowing any time gaps without compromising the string stability. On the contrary, as the ACC system has not C i-1 , the achievable stable time gap is longer, even close to the time gap of a manual driver. It explains why even with a high market penetration of ACC-equipped vehicles the traffic flow improvement effects are not visible.

K i G X i-1 - X i F i V c i-1 1/s h i - V i e i

B. Initial assumptions

The next assumptions have been considered in the analysis:

• Ego-vehicle (vehicle i) is at least the third vehicle in the string, meaning that there is a preceding one (vehicle i -1) whose V2V link is lost and a vehicle x (that could be the leader (x = 1), or some other if n > 3) broadcasting its information. • Maximum string length is limited to the speed difference value that makes drivers disengage ACC systems. Milanés et Al. [START_REF] Milanés | Cooperative adaptive cruise control in real traffic situations[END_REF] shows that this difference is equal to 5m/s and this value is reached in a six ACC-equipped vehicles string. V2V studies using Dedicated Short Range Communication (DSRC) guarantees up to 300m [START_REF] Armstrong | Dedicated short range communications (dsrc) home[END_REF], which also matches with the six-vehicles string assumption to guarantee communication even between leading and the last vehicle in the string. • The communicating vehicle ahead is in the same lane.

This can be known from a digital map and can be part of the transmitted information. • Communication link between vehicles is considered ideal.

No delay is present in the communication system.

III. YOULA-KUCERA CONTROLLER MODIFICATION

This section introduces the mathematical basis in which YK relies on to ensure stable interpolation between different controllers. The change between both controllers is made through the parameterization of all the controllers that stabilizes a given plant [START_REF] Kucera | A method to teach the parameterization of all stabilizing controllers[END_REF]. At the same time, the parameterization is based on the doubly coprime factorization [START_REF] Zhou | Robust and Optimal Control[END_REF].

The notation is standard; RH pxm ∞ is the real stable transfer function space which maps m-dimensional inputs into pdimensional outputs.

y = Gu, y ∈ R p and u ∈ R m ( 6
)
where G is the plant; y the measurement vector corresponding to the output and u the control signal. The dimension of both signals is

p = m = 1.
The plant could be stabilized by any appropriate LTI controller K i . Assuming that K 0 corresponds to CACC ST G and K 1 to CACC LT G , it will be possible to do a stable interpolation between them thanks to the YK parameterization of all controllers.

The plant G and both controllers K 0 and K 1 need to be factorized to apply YK theory.

A. Doubly coprime factorization

Double coprimeness refers to the idea of being right and left coprime [START_REF] Nett | A connection between state-space and doubly coprime fractional representations[END_REF].

Definition 1: Two different matrices M and N are right coprimes over RH pxm ∞ if they have the same number of columns and if matrices X r and Y r exist such that:

X r Y r M N = X r M +Y r N = I (7)
Definition 2: Two different matrices M and Ñ are left coprimes over RH pxm ∞ if they have the same number of rows and if matrices X l and Y l exist such that:

M Ñ X l Y l = MX l + ÑY l = I (8)
These coprime factorizations should be such that G and K i are:

G = NM -1 = M-1 Ñ K i = U i V -1 i = Ṽi -1 Ũi (9)
At the same time, these coprime factorizations

U i ∈ RH mxp ∞ , Ũi ∈ RH mxp ∞ , V i ∈ RH pxp ∞ , Ṽi ∈ RH pxp ∞ , N ∈ RH pxm ∞ , Ñ ∈ RH pxm ∞ , M ∈ RH pxm ∞ and M ∈ RH pxm ∞
satisfy the double Bézout's identity [START_REF] Pommaret | Generalized bezout identity[END_REF].

Ṽi -Ũi -Ñ M

M U i N V i = M U i N V i Ṽi -Ũi -Ñ M = I 0 0 I (10)
Coprime factorization related to any stabilizing controller in state space can be obtained if the theorem shown below is fulfilled [START_REF] Ishihara | Doubly coprime factorizations related to any stabilizing controllers in state space[END_REF]. Otherwise, the search of these factors could become complex.

Theorem 1: Consider a plant in state space representation as G(s) = C(sI -A) -1 B with A, B,C stabilizable and detectable and a stabilizing controller

K i (s) = C k,i (sI -A k,i ) -1 B k,i + D k,i . F and F k,i should be chosen such that A+BF and A k,i +B k,i F k,i ∈ RH pxm ∞ .
Then the matrices given by

M U i N V i =     A + BF 0 -B 0 0 A k,i + B k,i F k,i 0 B k,i -F C k,i + D k,i F k,i I D k,i -C F k,i 0 I     (11) Ṽi -Ũi -Ñ M =     A + BD k,i C BC k,i -B BD k,i B k,i C A k,i 0 B k,i F -D k,i C -C k,i I -D k,i C -F k,i 0 I     (12)
satisfy ( 9) and [START_REF] Arem | The impact of cooperative adaptive cruise control on traffic-flow characteristics[END_REF].

Notice how M and N remain the same for the plant G, as uncertainties are not taken into account in the model; while U 0 and V 0 will change to U 1 and V 1 when switching from the initial controller K 0 to the final one K 1 .

B. YK parameterization of all stabilizing controllers

YK parameterization provides all stabilizing controllers for a given plant G by interconnecting an initial controller with a parameter Q, called YK parameter, which can be any stable system with appropriate dimensions:

K = (U 0 + MγQ)(V 0 + NγQ) -1 = = ( Ṽ0 + γQ Ñ) -1 ( Ũ0 + γQ M), Q ∈ RH pxm ∞ (13)
From the general description of any stabilizing controller in equation [START_REF] Van Nunen | Robust model predictive cooperative adaptive cruise control subject to v2v impairments[END_REF], it is possible to do interpolation between controllers K 0 to K 1 online, without losing stability. To preserve a good transition between both controllers, the parameter Q should be regulated by a scalar factor γ ∈ [0, 1]. Mixed behaviours between both controllers can be obtained by applying different values of γ. But first, the value of Q should be obtained as follows (for further details [START_REF] Niemann | An architecture for implementation of multivariable controllers[END_REF]):

V 0 -1 U 0 Q N M - - γ K 0 r G K (γ)
u y Fig. 3. Control structure for modifying controllers online. Theorem 2: Let NM -1 be a coprime factorization of the model G and K 0 = U 0 V -1 0 = Ṽ0 -1 Ũ0 a stabilizing controller.

K i G X i-1 - X i F i /K i V c i-1 + 1/s G h i - K i ext V i
A second controller is given

K 1 = U 1 V -1 1 = Ṽ1 -1 Ũ1 . If G * K 1 ∈ RH ∞ ∧ Ṽ1 Ṽ0 -1 V 1 ∈ RH ∞ ( 14
)
It is equivalent to the existence of Q, so all the stabilizing controllers are described by:

K 1 = (V 0 + Q Ñ) -1 (U 0 + Q M) ( 15
)
where Q can be calculated as:

Q = Ũ1 -Ṽ1 Ṽ0 -1 Ũ0 (16)
This theorem allows obtaining a control structure able to do stable interpolation between different controllers. Theorem 2 is used for getting the control structure shown in Fig. 3.

When doing controller interpolation, γ plays a key role. It regulates the different levels of activation of the YK parameter Q. γ may vary from 0 to 1, being 0 a 100% contribution of K 0 and 1 a 100% contribution of K 1 .

The stability in the transition is ensured when Q ∈ RH ∞ ; so any linear combination of Q and γ will provide stable responses [START_REF] Niemann | An architecture for implementation of multivariable controllers[END_REF].

K i ext k p,i k d,i h i CACC ST G K 0 ext k p,0 k d,0 h 0 = h CACC = 0.6s CACC LT G K 1 ext k p,1 k d,1 h 1 = h ACC = 1.5s
IV. ACACC ACACC is composed by two different CACC controllers: CACC ST G and CACC LT G . Time gap values h CACC and h ACC are chosen according to the general public accepted time gaps for CACC and ACC systems [START_REF] Nowakowski | Cooperative adaptive cruise control: Driver acceptance of following gap settings less than one second[END_REF]. For CACC, the shortest gap is set at 0.6s, while for ACC is set to 1.5s.

For each of these controllers, the PD-based CACC controller seen in Fig. 2 is reformulated in Sec. IV.A to include the time gap h i into the controller. This modification of the classical structure allows stable interpolation between controllers with different h i , k p,i and k d,i . By taking the general theory in Sec. III a control structure for making stable transitions between CACC ST G and CACC LT G is shown in Sec. IV.B.

Once the control structure is obtained, a different percentage of each of the controllers is applied depending on the traffic situation. The percentage is chosen through the YK gamma γ.

The algorithm for discriminating if the information received via V2V is coming either from vehicle i -1 or from vehicle x in the string is out of the scope of this paper (see [START_REF] Bu | Design and field testing of a cooperative adaptive cruise control system[END_REF] for a possible implementation).

In Sec. IV.C, γ is designed to be modified depending on the communicated velocity V x and the velocity of the preceding vehicle received through the onboard sensor system V s i-1 .

A. Modification of classical CACC structure

The classical PD-based CACC controller is modified in Fig. 4 to include the time gap h i into the controller. In that way, it will be possible to change the time gap of the system when doing controllers interpolation. The extended controller K i ext results in the transfer function [START_REF] Shahruz | Design of controllers for linear parametervarying systems by the gain scheduling technique[END_REF]. The feedforward filter is changed to F i /K i .

CACC ST G has the index i = 0, while CACC LT G has the index i = 1. The information related to both controllers is summarized in Table I. Notice how the values of k p,0 , k p,1 , k d,0 and k d,1 will depend on the accepted time gaps by drivers and the vehicle model G.

K i ext (s) = K i 1 + Gh i K i = k p,i + k d,i s 1 + Gh i (k p,i + k d,i s) (17) V 0 -1 U 0 Q N M - - γ K 0 ext X i-1 G F i /K i 1/s V c 1 X i K ext V i (γ)

B. Control structure for CACC controller modification

Once CACC ST G and CACC LT G are defined, is possible to to have hybrid behaviours between both, without losing stability by obtaining the parameter Q with Theorem 2. Figure 5 modifies the control structure in Fig. 3 to allow CACC controller modification. Please notice how the structure only differs in adding the filtered communication link with the closest V2Vequipped vehicle, in that case, communication link C x .

C. YK gamma tuning

Once the YK control structure for CACC is obtained, a different percentage of each of the controllers is applied depending on the traffic situation. γ = 0 when C i-1 exists, while a simple decision-making system is used when C i-1 does not exist, but there is communication with a vehicle x, C x . Gamma tuning is chosen in order to improve traffic flow:

γ = 0.033(V s -V x ) + 0.5 if | V x -V s i-1 |< 5m/s&∃C x & C i-1 0 if ∃C i-1
where V s i-1 is the velocity of the vehicle i -1 obtained through the on-board sensor systems, and V x is the velocity of vehicle x received through C x . A maximum speed oscillation of 5m/s is considered as the operation range for the present application.

Different traffic situations are pictured for understanding the tuning of γ when C i-1 . When vehicle i -1 has a behaviour similar to vehicle x, γ is modified to response faster to speed changes:

• If vehicle x accelerates, V x will be higher than V s i-1 . So, γ decreases, making the gap time shorter.

• If vehicle x brakes, V x will be lower than V s i-1 . So, γ increases, making the gap time longer.

• When both vehicle x and vehicle i -1 have similar velocities, γ results in an intermediate value between both controllers. The feedforward filter F is modified accordingly to the gamma range as:

F =      (0.033(V s i-1 -V x )) + 1.0) • F 0 /K 0 if | V x -V s i-1 |< 5m/s&∃C x & C i-1 F 0 /K 0 if ∃C i-1 V. SIMULATION RESULTS
This section presents the ACACC performance at high speeds (i.e. highway scenario). For simulation purposes, the vehicle model introduced in [START_REF] Milanés | Cooperative adaptive cruise control in real traffic situations[END_REF] is used. It is a second-order response with a time delay coming from a Nissan Infiniti M56. For the sake of clarity, vehicle transfer function is shown in [START_REF] Stilwell | Interpolation of observer state feedback controllers for gain scheduling[END_REF].

G(s) = 1.136 s 2 + 1.067s + 1.1385 e -0.287s (18)
The algorithm performance is tested by using as vehicle i-1 a non V2V-equipped vehicle. It is modelled as a human driverusing the Intelligent Driver Model (IDM) [START_REF] Kesting | Enhanced intelligent driver model to access the impact of driving strategies on traffic capacity[END_REF]-or equipped with an ACC system.

A. Control system stability

CACC ST G and CACC LT G controller parameters are shown in Table II. Control system stability for the designed ACACC is studied in function of γ. A necessary and sufficent condition to ensure control system stability for a feedback/feedforward control loop as the one in Fig. 5 is:

I -K ext (γ) -G/s I -1 ∈ RH ∞ = = K ext (γ)G/s 1 + K ext (γ)G/s = K 0 ext +M 0 γQ 1+N 0 γQ G/s 1 + K 0 ext +M 0 γQ 1+N 0 γQ G/s ∈ RH ∞ with F 0 /K 0 ∈ RH ∞ (19)
As the feedforward filter is already stable, poles of the resulting transfer function in equation ( 19) are analysed for different values of γ. Closed-loop (CL) poles are in Table III. CL poles during the transition are the combination of CL poles of (G/s, K 0 ext ) and (G/s, K 1 ext ). Control system stability is ensured for every value of γ.

String stability is also analysed. First, when there is communication with the preceding vehicle (see equation ( 3)): Fig. 6 represents the stability surface for the designed ACACC in function of the γ value. One can appreciate how the infinity norm X i /X i-1 is always below the unity for any value of γ, proving string stability. Second, when there is no communication with the preceding vehicle, but another further on: Control velocity V c

x (with x = 1 in the case of a string of three vehicles) is received instead of V c i-1 . The designed YK gamma tuning is employed in this situation. Figure 7 considers difference between V s i-1 and V 1 up to 15m/s. As ACACC application range is up to 5m/s, string stability is ensured even if communication with the preceding vehicles is not available.

K i ext k p,i k d,i h i CACC ST G K 0 ext k p,0 = 0.45 k d,0 = 0.25 h 0 = h CACC = 0.6s CACC LT G K 1 ext k p,1 = 0.45 k d,1 = 0.25 h 1 = h ACC = 1.5s

B. Using IDM as preceding non V2V-equipped vehicle

First simulations consider a human-driven vehicle (modelled by the IDM) as vehicle i-1 without communication. IDM is a well-known car-following model in the traffic flow simulation literature [START_REF] Kesting | Enhanced intelligent driver model to access the impact of driving strategies on traffic capacity[END_REF]. It defines a target acceleration depending on desired velocity, desired time gap and comfortable deceleration. The parameters in [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] are used as reference values.

Figure 8 depicts the performance of the ACACC controller in comparison with ACC and CACC (assuming V2V capabilities in the IDM vehicle, only in the CACC case) controllers The value of γ is always around 0.5, changing with acceleration or braking changes. The performance of the ACACC system is a hybrid response between CACC ST G and CACC LT G . Its response comes earlier to changes in leader speed V x than the ACC/CACC systems. For instance, when vehicle 1 is braking at 300s, the ACACC system brakes practically at the same time, while ACC and CACC react 5s later. Besides, the system is more damped while speed changes take place, better preserving string stability. Finally, the main objective of the present paper is also fulfilled, instead of degrading to an ACC system, ACACC takes C 1 for making the string of vehicles tighter, improving highway capacity. The relative distance corresponds to a spacing policy between CACC ST G and CACC LT G responses. Figure 9 shows a second simulation using IDM. The desired velocity parameter from [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] is change from 33.33m/s to 28m/s, emulating a very slow vehicle dynamic on the non V2V-equipped vehicle i-1 when tracking changes in V x . Since the IDM model is not correctly following vehicle 1 response, the ACACC system has a closer behaviour to CACC LT G . This behaviour is appreciated on the lower plot of Fig. 9.

C. Using an ACC-equipped vehicle as preceding non V2Vequipped one

Milanés and Shladover [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] compare the car following performance among ACC and CACC systems for a Nissan Infiniti M56s. Experimental results are used to obtain carfollowing models for representing the production ACC and the new designed CACC controller, highlighting how the production ACC results string unstable, amplifying the speed changes of the preceding vehicle. A simulation using the ACC/CACC car-following models presented in [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] is carried out. Specifically, the more the ACC-equipped vehicles in the string, the bigger the amplification. Here, instead of degrading the system to another ACC, the system is changed to ACACC with C 1 .

Figure 10 depicts the performance of the ACACC controller in comparison with ACC and CACC controllers when an ACC-equipped vehicle is the vehicle i -1. Vehicle i -1 is a V2Vequipped vehicle that loses its communication capability at 220s. A string of three vehicles is considered. For notation, vehicle 1 (solid magenta line) is a CACC-equipped vehicle following a speed reference; vehicle i -1 (solid cyan line) is the one that starts the platoon in the second position with an ACC as controller; and, vehicle i or follower is either an ACCequipped vehicle (solid red line), a CACC-equipped vehicle (solid green line, considering that vehicle i -1 is always V2Vequipped) or an ACACC-equipped vehicle (solid blue line).

ACACC works as a classical CACC system while the communication with the preceding vehicle is available. Once the communication is lost, ACACC system significantly reduces speed oscillation introduced by the ACC-equipped vehicle, providing string stability but also increasing traffic flow by reducing inter-vehicle distances.

Finally, limits of the ACACC applications are tested. As already said, 5m/s is used as the maximum acceptable range for using ACACC. Specifically, Fig. 11 depicts the performance of the ACACC controller in comparison with ACC and CACC controllers when four ACC-equipped vehicles are between vehicle 1 and vehicle i. So, a string of six vehicles is considered. Speed difference between vehicle 1 and vehicle i-1 reach almost 5m/s. V2V communication is only available between vehicle 1 and i. For notation, vehicle 1 (solid magenta line) is a CACC-equipped vehicle following a speed reference; vehicle i -1 (solid cyan line) is the one that starts the platoon in the fifth position with an ACC as controller; and, vehicle i or follower is either an ACC-equipped vehicle (solid red line), a CACC-equipped vehicle (solid green line, considering that vehicle i-1 is V2V-equipped) or an ACACC-equipped vehicle (solid blue line).

Even with 4 ACC-equipped vehicles without communications between leader and ego vehicles, ACACC preserves its ability to reduce speed oscillations and intervehicle-distance, improving consequently the traffic flow.

These simulations represent the closer behaviour to real traffic environment with a high penetration of ACC-equipped and CACC-equipped vehicles, giving an insight into the potential benefits of ACACC.

VI. EXPERIMENTAL RESULTS

Three cycabs were used as experimental vehicles. These vehicles are mobile platforms used in several research labs conceived for urban applications. From a kinematic point of view, cycabs can turn their rear wheels as a linear function of the steering angle of the front wheels [START_REF] Sekhavat | The cycab robot: a differentially flat system[END_REF]. The vehicle is controlled by a commanded velocity. The velocity is limited up to 4m/s , because it is mainly designed for crowded areas where higher speeds can lead to unsafe situations. As ACACC has been designed for high-velocities, the proposed system is scaled to these low-speed platforms for having an experimental proof-of-the-concept.

The vehicle dynamic model for evaluating the performance of ACACC is presented in equation [START_REF] Niemann | An architecture for implementation of multivariable controllers[END_REF]. Cycab real velocity responses have a characteristic overshoot preserved in the model.

G(s) = 1 0.8768s 2 + 1.252s + 1 (20)
A. Control system stability γ = 0.1 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.2 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.3 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.4 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.5 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.6 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.7 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.8 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 0.9 [-6.28e 7 , -6.28e 7 , 0, 0, -0.78 ± 1.24, -0.88 ± 1.70i] γ = 1.0 [-6.28e 7 , 0, -0.88 ± 1.70i] Frequency (rad/s) Fig. 13. String stability magnitude in dB for a cycab when communication with preceding vehicle is not available, but it is with the leader vehicle. |

K i ext k p,i k d,i h i CACC ST G K 0 ext k p,0 = 1.5 k d,0 = 0.2 h 0 = h CACC = 0.6s CACC LT G K 1 ext k p,1 = 1.5 k d,1 = 0.2 h 1 = h ACC = 1.5s
V 1 -V s i-1 |= 1.5m/s..
Control system stability for the designed ACACC is studied in function of γ through equation [START_REF]Stability preserving interpolation methods for the synthesis of gain scheduled controllers[END_REF]. Controller gains in Table IV and model in equation [START_REF] Niemann | An architecture for implementation of multivariable controllers[END_REF] are considered in order to obtain the CL poles shown in Table V. Control system stability is preserved for every value of γ.

String stability is analysed. When there is communication with the preceding vehicle, Fig. 12 represents the stability surface for ACACC in function of the γ value. Its infinity norm is always below the unity for any value of γ, proving string stability. On the other hand, when there is no communication with the preceding vehicle but with a leader in a string of three vehicles, the bode's diagram in Fig. 13 considers a difference between V s i-1 and V 1 of 1.5m/s (a scale of 10 is considered in the YK gamma tuning for obtaining a low-velocity behaviour). String stability is again ensured even if communication with the preceding vehicle is not available.

B. Using IDM as preceding non V2V-equipped vehicle

The non V2V-equipped vehicle i -1 is a cycab with an IDM controller with the desired velocity of 33.33m/s. As IDM controllers depend on the ego-velocity and the velocity of its preceding vehicle, both signals are scaled with a factor of 10 for obtaining a low-velocity behaviour. In the same way, the velocities used in the tuning of γ are scaled for proving the correct performance of the approach.

Figure 14 depicts the performance of the scaled ACACC controller in comparison with ACC and CACC controllers when an IDM-equipped cycab is the vehicle i -1. The top graph plots the velocity responses of the three cycabs in the string. The second depicts the different accelerations of these vehicles during the test. The third graph shows the relative distances to the preceding vehicle. And the bottom graph represents how γ is changed over time, together with the real time gap.

For notation, vehicle 1 (solid magenta line) is a cycab with a speed profile that goes from 1m/s to 3m/s or vice-versa with increasing steps of acceleration/deceleration; vehicle i-1 (solid cyan line) is the second cycab in the string with an IDM controller; and, vehicle i or follower in third position is either ACC (solid red line), CACC (solid green line) or ACACC (solid blue line).

The performance of the ACACC is well-proved. The hybrid behaviour between CACC ST G and CACC LT G responses of the ACACC controller is preserved. Its response comes earlier when vehicle 1 is accelerating or braking. And the relative distance to vehicle i -1 is shorter than with the classical degradation to ACC, improving the capacity of the road.

Nevertheless, it is noteworthy to mention how at lowvelocities the adapted time gap by the system is more sensitive to any velocity change, as well to the standstill distance. It explains why the variation of the real time gap time through γ becomes more complicated for the low-velocity test (see the bottom part in Fig. 14).

C. Using an ACC-equipped vehicle as preceding non V2Vequipped one

A scenario where vehicle i -1 is an ACC-equipped cycab is carried out. Please notice that the ACC shown here is string stable. The performance of the ACACC controller can be accomplished with a stable ACC without endangering the experimental platforms.

Figure 15 depicts the performance of the ACACC controller compared to ACC and CACC when the non V2V-equipped vehicle is an ACC-equipped cycab. For notation, vehicle 1 (solid magenta line) is a cycab with the speed profile seen in the previous section; vehicle i -1 (solid cyan line) is the second cycab with an ACC controller; and vehicle i or follower in the third position is either ACC (solid red line), CACC (solid green line) or ACACC (solid blue line).

These low-speed tests demonstrate the good performance of the proposed system. ACACC controller responses faster to speed changes, while preserving better the string stability and improving the highway capacity with a shorter distance to vehicle i -1.

VII. CONCLUSIONS

This paper explores the benefit of using V2V-equipped vehicle information from a vehicle ahead in the string when the preceding vehicle is a non-equipped one to improve traffic flow. The proposed system provides string stable responses and, more interestingly, reduce inter-vehicle distances, providing tighter strings. System performance is evaluated using either a human driven vehicle (emulated by the IDM) or an ACC-equipped vehicle as preceding one.

The stable interpolation between CACC controllers with different time gaps is ensured thanks to the use of the YK parameterization. The string stability of the control structure is also demonstrated.

A comparative analysis using [START_REF] Milanés | Modeling cooperative and autonomous adaptive cruise control dynamic responses using experimental data[END_REF] models for ACC and CACC vehicles at high speeds and the IDM model as human driver vehicle is carried out. ACACC system exhibits a good performance, improving ACC response by providing a tighter string of vehicles, a faster response to speed changes and a better preservation of the string stability. These results are scaled to the low-velocity experimental platform for having an experimental proof-of-concept. The performance is wellproved, improving relative distance to the preceding vehicle and obtaining faster responses. Highway capacity improvement when using ACACC will be quantified and studied in function of the vehicles market penetration, to have a more fair comparison with results in [START_REF] Shladover | Impacts of cooperative adaptive cruise control on freeway traffic flow[END_REF] and [START_REF] Werf | Effects of adaptive cruise control systems on highway traffic flow capacity[END_REF].

Non-ideal V2V communication link needs to be studied in a future, analyzing the relationship between stability and delay. This would need another controller configuration that could be included into the YK parameterization. Additionally, when considering different braking capabilities in each of the vehicles in the string, a collision avoidance system is needed, which represents a different control system to be developed and included in the YK approach.

Fig. 1 .

 1 Fig. 1. ACACC scenario.

Fig. 2 .

 2 Fig. 2. Classical CACC control structure.

Fig. 4 .

 4 Fig. 4. Modified CACC control structure.

Fig. 5 .

 5 Fig. 5. YK control structure for modifying CACC controllers online.

γFig. 6 .

 6 Fig. 6. String stability's surface for a Nissan Infiniti M56s depending on YK γ. Preceding vehicle's communication is available.

Fig. 7 .

 7 Fig. 7. String stability magnitude in dB for a Nissan Infiniti M56 when communication with preceding vehicle is not available, but it is with another vehicle further on the string. | V x -V s i-1 |= 15m/s.

Fig. 8 .

 8 Fig. 8. Simulation results comparison among car-following policies using CACC, ACC and ACACC controllers when a non V2V-equipped vehicle is in front. The non V2V-equipped is a conventional driver emulated by a IDM with a desired velocity of 33.33m/s.

Fig. 9 .

 9 Fig.9. Simulation results comparison among car-following policies using CACC, ACC and ACACC controllers when a non V2V-equipped vehicle is in front. The non V2V-equipped is a conventional driver emulated by an IDM with a desired velocity of 28m/s.

Fig. 10 .Fig. 11 .

 1011 Fig. 10. Simulation results comparison among car-following policies using CACC, ACC and ACACC controllers when communication is lost with a V2V-equipped vehicle. Preceding vehicle is an ACC-equipped vehicle with communication capabilities. Communication is lost at 220s.

Fig. 12 .

 12 Fig. 12. String stability's surface for a Cycab depending on YK γ.

  s Vi-1=1.5 m/s V1=1.5 m/s Vi-1=3 m/s Magnitude Xi/Xi-1 (dB)

Fig. 14 .Fig. 15 .

 1415 Fig.[START_REF] Savaia | A receding horizon switching control resilient to communication failures for connected vehicles[END_REF]. Experimental results comparison among car-following policies using CACC, ACC and ACACC controllers when a non V2V-equipped cycab is in front. The non V2Vequipped is a cycab with an IDM with a desired velocity of 33.33m/s.