D Dragičević
email: d.dragicevic@unsw.edu.au.

G Froyland
email: g.froyland@unsw.edu.au.

Cecilia González-Tokman
email: cecilia.gt@uq.edu.au.

S Vaienti
email: vaienti@cpt.univ-mrs.fr.

Almost sure invariance

published or not. The documents may come

Introduction

The objective of this note is to prove the almost sure invariance principle (ASIP) for a large class of random dynamical systems. The random dynamics is driven by an invertiblemeasure preserving transformation σ of (Ω, F, P) called the base transformation. Trajectories in the phase space X are formed by concatenations

f n ω := f σ n-1 ω • • • • • f σω • f ω
of maps from a family of maps f ω : X → X, ω ∈ Ω. For a systematic treatment of these systems we refer to [START_REF] Arnold | Random Dynamical Systems[END_REF]. For sufficiently regular bounded observables ψ ω : X → R, ω ∈ Ω, an almost sure invariance principle guarantees that the random variables ψ σ n ω • f n ω can be matched with trajectories of a Brownian motion, with the error negligible compared to the length of the trajectory. In the present paper, we consider observables defined on some measure space (X, m) which is endowed with a notion of variation. In particular, we consider examples where the observables are functions of bounded variation or quasi-Hölder functions on a compact subset X of R n . We emphasize that our setting is quite similar to that in [START_REF] Buzzi | Exponential Decay of Correlations for Random Lasota-Yorke Maps[END_REF], where the maps f ω are called random Lasota-Yorke maps.

In a more general setting and under suitable assumptions, Kifer proved in [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF] central limit theorems (CLT) and laws of iterated logarithm; we will briefly compare Kifer's assumptions with ours in Remark 2 below. In [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF]Remark 2.7], Kifer claimed without proof (see [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF]Remark 4.1]) a random functional CLT, i.e. the weak invariance principle (WIP), and also a strong version of the WIP with almost sure convergence, namely the almost sure invariance principle (ASIP), referring to techniques of Philip and Stout [START_REF] Philip | Almost sure invariance principle for sums of weakly dependent random variables[END_REF].

Here we present a proof of the ASIP for our class of random transformations, following a method recently proposed by Cuny and Merlèvede [START_REF] Cuny | Strong invariance principles with rate for "reverse" martingales and applications[END_REF]. This method is particularly powerful when applied to non-stationary dynamical systems; it was successfully used in [START_REF] Haydn | Almost sure invariance principle for sequential and non-stationary dynamical systems[END_REF] for a large class of sequential systems with some expanding features and for which only the CLT was previously known [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF]. We stress that ω-fibered random dynamical systems discussed above are also non-stationary since we use ω-dependent sample measures (see below) on the underlying probability space.

The technique of Cuny and Merlèvede is based on martingale approximation; it was shown in [START_REF] Haydn | Almost sure invariance principle for sequential and non-stationary dynamical systems[END_REF] how to satisfy one of the main assumptions in [START_REF] Cuny | Strong invariance principles with rate for "reverse" martingales and applications[END_REF] by using a result of Sprindzuk [START_REF] Sprindzuk | Metric theory of Diophantine approximations[END_REF], which basically consists of getting an almost sure bound when the latter is known to hold in the L 1 norm. To prove such a result we also need two other ingredients: (i) the error in the martingale approximation must be bounded in a suitable Banach space; (ii) the quenched correlations with respect to the sample measures must decay with a summable rate.

We now compare our assumptions and results with those in Kifer's paper [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF]. Kifer used a martingale approximation, but the martingale approximation error in [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF] is given in terms of an infinite series (see the error g ω in equation (4.18) in [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF]), which appears difficult to estimate under general assumptions. Instead, our martingale approximation error term is explicitly given in terms of a finite sum (see ( 14)), and as mentioned above we can bound it easily. Furthermore, Kifer invoked a rate of mixing, but to deal with it he assumed strong conditions (φ-mixing and α-mixing), which are very hard to check on concrete examples. We use instead quenched decay of correlations on a space of regular observables, for example, bounded variation or quasi-Hölder and L ∞ functions (exponential decay was shown by Buzzi [START_REF] Buzzi | Exponential Decay of Correlations for Random Lasota-Yorke Maps[END_REF]), with an addition: the constant that scales the norm of the observable in the decay rate is independent of the noise ω; we can then satisfy the hypotheses of Sprindzuk's result. Further comparisons are deferred to Remark 2.

The rate which we obtain by approximating our process with a sum of i.i.d Gaussian variables (the content of the ASIP) is of order n 1/4 , which is, up to logarithmic corrections, a rate previously obtained for deterministic uniformly expanding systems [START_REF] Gouëzel | Almost sure invariance principle for dynamical systems by spectral methods[END_REF].

Preliminaries and statement of the main results

2.1.

Preliminaries. We introduce in this section the fiber maps and the associated function spaces which we will use to form the random concatenations. We will call them random expanding transformations, or random Lasota-Yorke maps. We will refer to and use the general assumptions for these maps as proposed by Buzzi [3] in order to use his results on quenched decay of correlations. However, we will strengthen a few of those assumptions with the aim of obtaining limit theorems. Our additional conditions are similar to those called Dec and Min in the paper [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF], and which were used to establish and recover a property akin to quasi-compactness for the composition of transfer operators.

Let (Ω, F, P) be a probability space and let σ : Ω → Ω be an invertible P-preserving transformation. We will assume that P is ergodic. Moreover, let (X, B) be a measurable space endowed with a probability measure m and a notion of a variation var :

L 1 (X, m) → [0, ∞] which satisfies the following conditions: (V1) var(th) = |t| var(h); (V2) var(g + h) ≤ var(g) + var(h); (V3) h ∞ ≤ C var ( h 1 + var(h)) for some constant 1 ≤ C var < ∞; (V4) for any C > 0, the set {h : X → R : h 1 + var(h) ≤ C} is L 1 (m)-compact; (V5) var(1 X ) < ∞, where 1 X denotes the function equal to 1 on X; (V6) {h : X → R + : h 1 = 1 and var(h) < ∞} is L 1 (m)-dense in {h : X → R + : h 1 =
1}. (V7) there exists K var < ∞ such that var(f g) + f g 1 ≤ K var (var(f ) + f 1 )(var(g) + g 1 ), for every f, g ∈ BV .

(

) (V8) for any f ∈ L 1 (X, m) such that essinf f > 0, we have var(1/f ) ≤ var(f ) (essinf f ) 2 . 1
We denote by BV = BV (X, m) the space of all h ∈ L 1 (X, m) such that var(h) < ∞. It is well known that BV is a Banach space with respect to the norm

h BV = var(h) + h 1 .
On several occasions we will also consider the following norm

h var = var(h) + h ∞ , on BV which (although different) is equivalent to • BV . Let f ω : X → X, ω
∈ Ω be a collection of mappings on X. The associated skew product transformation τ : Ω × X → Ω × X is defined by

τ (ω, x) = (σ(ω), f ω (x)).
(

) 2
Each transformation f ω induces the corresponding transfer operator L ω acting on L 1 (X, m) and defined by the following duality relation

X (L ω φ)ψ dm = X φ(ψ • f ω ) dm, φ ∈ L 1 (X, m), ψ ∈ L ∞ (X, m).
For each n ∈ N and ω ∈ Ω, set

f n ω = f σ n-1 (ω) • • • • • f ω and L n ω = L σ n-1 (ω) • • • • • L ω .
We say that the family of maps f ω , ω ∈ Ω (or the associated family of transfer operators

L ω , ω ∈ Ω) is uniformly good if: (H1) The map (ω, x) → (L ω H(ω, •))(x) is P × m-measurable for every P × m-measurable function H such that H(ω, •) ∈ L 1 (m) for a.e. ω ∈ Ω; (H2) There exists C > 0 such that L ω φ BV ≤ C φ BV
for φ ∈ BV and P a.e. ω ∈ Ω. (H3) For P a.e. ω ∈ Ω, sup

n≥0 φ n+1 • f σ n (ω) BV < ∞,
whenever {(φ n ) n≥0 } ⊂ BV and sup n φ n BV < ∞. (H4) There exist K, λ > 0 such that L n ω φ BV ≤ Ke -λn φ BV , for n ≥ 0, P a.e. ω ∈ Ω and φ ∈ BV such that φ dm = 0. (H5) There exists c > 0 such that L n ω 1 X ≥ c, for P a.e. ω ∈ Ω and n ∈ N.

Using (H1), (H2), and (H4) we can prove the existence of a unique random ACIM h : Ω × X → R, with uniformly bounded fibres h ω .

Proposition 1. Let f ω , ω ∈ Ω be a uniformly good family of maps on X. Then there exist a unique measurable and nonnegative function h : Ω × X → R with the property that

h ω := h(ω, •) ∈ BV , h ω dm = 1, L(h ω ) = h σ(ω) for a.e. ω ∈ Ω esssup ω∈Ω h ω BV < ∞. (3)
Proof. Let

Y = v : Ω × X → R : v measurable, v ω := v(ω, •) ∈ BV and esssup ω∈Ω v ω BV < ∞ .
Then, Y is a Banach space with respect to the norm

v ∞ := esssup ω∈Ω v ω BV .
Moreover, let Y 1 be the set of all v ∈ Y such that v ω dm = 1 and v ω ≥ 0 for a.e. ω ∈ Ω.

It is easy to verify that Y 1 is a closed subset of Y and thus a complete metric space. We define a map L :

Y 1 → Y 1 by (L(v)) ω = L σ -1 (ω) v σ -1 (ω) , ω ∈ Ω, v ∈ Y 1 .
Note that it follows from (H2) that

L(v) ∞ = esssup ω∈Ω (L(v)) ω BV ≤ C esssup ω∈Ω v σ -1 (ω) BV = C v ∞ .
Furthermore,

(L(v)) ω dm = L σ -1 (ω) v σ -1 (ω) dm = v σ -1 (ω) dm = 1,
for a.e. ω ∈ Ω. Hence, L is well-defined. Similarly,

L(v) -L(w) ∞ ≤ C v -w ∞ , for v, w ∈ Y 1
which shows that L is continuous. Choose n 0 ∈ N such that Ke -λn 0 < 1. Take arbitrary v, w ∈ Y 1 and note that by (H4),

L n 0 (v) -L n 0 (w) ∞ = esssup ω∈Ω L n 0 σ -n 0 (ω) (v σ -n 0 (ω) -w σ -n 0 (ω)) BV ≤ Ke -λn 0 esssup ω∈Ω v σ -n 0 (ω) -w σ -n 0 (ω) BV = Ke -λn 0 v -w ∞ .
Hence, L n 0 is a contraction on Y 1 and thus has a unique fixed point h ∈ Y 1 . Set

h ω := 1 n 0 hω + 1 n 0 L σ -1 (ω) ( hσ -1 (ω) ) + . . . + 1 n 0 L n 0 -1 σ -(n 0 -1) ω ( hσ -(n 0 -1) ω ), ω ∈ Ω.
Then, h is measurable, nonnegative, h ω dm = 1 and a simple computation yields L(h ω ) = h σ(ω) . Finally, by (H2) we have that

esssup ω∈Ω h ω BV ≤ C n 0 -1 n 0 (C -1) esssup ω∈Ω hω BV < ∞.
Thus, we have established existence of h. The uniqueness is obvious since each h satisfying the assertion of the theorem is a fixed point of L and thus also of L n 0 which implies that it must be unique.

We note that [START_REF] Buzzi | measures for random Lasota-Yorke maps[END_REF][START_REF] Buzzi | Exponential Decay of Correlations for Random Lasota-Yorke Maps[END_REF] prove the above existence result with weaker control on the properties of f ω , and obtain existence of a random ACIM {h ω } ω∈Ω under less restrictive conditions. Indeed, those results don't require (H4) and in addition (H2) is allowed to hold with C = C(ω) such that log C ∈ L 1 (P).

We now describe a large class of examples of good families of maps f ω , ω ∈ Ω. We first show that they satisfy properties (H1)-(H3); this will crucially depend on the choice of the function space. We then give additional requirements in order for those maps to satisfy condition (H4), also called Dec in [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF] when applied to sequential systems, and condition (H5), named Min in [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF].

2.2. Examples. 2.2.1. Random Lasota-Yorke maps. Take X = [0, 1], a Borel σ-algebra B on [0, 1] and the Lebesgue measure m on [0, 1]. Furthermore, let var(g) = inf h=g(mod m) sup 0=s 0 <s 1 <...<sn=1 n k=1 |h(s k ) -h(s k-1 )|.
It is well known that var satisfies properties (V1)-(V8) with C var , K var = 1. For a piecewise

C 2 f : [0, 1] → [0, 1], set δ(f ) = essinf x∈[0,1] |f |. Consider now a finitely-valued, measurable map ω → f ω , ω ∈ Ω of piecewise C 2 maps on [0, 1] satisfying (H1) such that sup ω∈Ω N (f ω ) =: N < ∞, inf ω∈Ω δ(f ω ) =: δ > 1, and sup ω∈Ω |f ω | ∞ =: D < ∞.
It is proved in [START_REF] Buzzi | Exponential Decay of Correlations for Random Lasota-Yorke Maps[END_REF] that the family f ω , ω ∈ Ω satisfies (H2) with

C = 4 N δ ∨ 1 D δ 2 ∨ 1 1 δ ∨ 1 , (4)
where for any two real-valued functions g 1 and g 2 , g 1 ∨ g 2 = max{g 1 , g 2 }, and (V8) has been used for the bound var(1/f ) ≤ D δ 2 . We note that since N < ∞, condition (H3) holds.

The uniform decay rate (H4) has been treated in Propositions 2.10 and 2.11 in [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF]. There, sufficient conditions were stated for sequential dynamical systems, but these conditions can be can be easily adapted to our random setting. Conze and Raugi [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF] propose two types of conditions, either of which yields (H4). The first local type requires the existence of a fiber map, say f 0 , whose transfer operator L 0 is quasi-compact and exact. Defining a distance between two transfer operators [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF] that there is a neighborhood U 0 of L 0 such that all allowed concatenations of transfer operators drawn from U 0 verify (H4) as well.

L 1 , L 2 , as d(L 1 , L 2 ) := sup φ∈BV, φ BV ≤1 ||L 1 φ -L 2 φ|| 1 , it was proved in
To establish (H4) for the second, more general, nonlocal type of random dynamics we require two conditions:

• Random covering: Let A ω denote the collection of intervals of monotonicity for the map f ω and define A n ω = n-1 j=0 (f j ω ) -1 A σ j ω . We will say that the random Lasota-Yorke maps {f ω } ω∈Ω are covering if for each n ≥ 0, ω ∈ Ω, and J ∈ A n ω there is an

n 0 such that f n 0 ω (J) = [0, 1]. • Uniform Lasota-Yorke inequality: There exist n ∈ N, 0 < ρ < 1 and B > 0 such that for a.e. ω ∈ Ω, L n ω f BV ≤ ρ f BV + B f 1 .
For a fixed covering Lasota-Yorke map, Liverani [12, Theorem 3.6] established exponential decay of correlations for observables of bounded variation, using cone techniques and the property that the unique random invariant density is uniformly bounded below [START_REF] Liverani | Decay of correlations for piecewise expanding maps[END_REF]Lemma 4.2]. His results, in particular [12, Lemma 3.5], which determines the rate of correlation decay, are directly applicable in our setting of random composition of finitely many Lasota-Yorke maps. This ensures exactness of the sequences (f σ j ω ) j≥0 for every ω ∈ Ω. That is, for every φ ∈ BV with φ dm = 0, lim n→∞ L n ω φ 1 = 0; see for example the proof of Proposition 3.6 in [START_REF] Froyland | Stability and approximation of random invariant densities for Lasota-Yorke map cocycles[END_REF].

The proof of (H4) now follows as the proof of [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF]Proposition 2.11]. Indeed exactness, together with the compactness condition (for every sequence (L ω j ) j∈N , there exists a convergent subsequence) ensures that for every ε 0 > 0 and every ω ∈ Ω there exists q ∈ N such that for every φ ∈ BV with φ dm = 0 and every j ∈ N, L q σ j ω φ 1 ≤ ε 0 φ BV , via a diagonal argument. This property, combined with the uniform Lasota-Yorke inequality implies (H4) essentially as in the proof of [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF]Proposition 2.7].

With the random covering and uniform Lasota-Yorke conditions introduced above, and assuming esssup ω∈Ω |f ω | ∞ ≤ C , it follows from Proposition 2 in [START_REF] Aimino | Concentration inequalities for sequential dynamical systems on the unit interval, Ergodic Theory Dynam[END_REF] that (H5) holds. We note that [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF] demonstrated that (H5) holds for compositions of β-transformations with β selected from an appropriate interval of values, and that [START_REF] Haydn | Almost sure invariance principle for sequential and non-stationary dynamical systems[END_REF] stated similar sufficient conditions for (H5) for sequences of Lasota-Yorke maps that are translations of a fixed Lasota-Yorke map or small perturbations of a fixed Lasota-Yorke map.

2.2.2. Random piecewise expanding maps. In higher dimensions, the properties (V1)-(V8) can be checked for a so-called quasi-Hölder space, which in particular is injected in L ∞ (condition (V3)) and has the algebra property (V7). Originally developed by Keller [START_REF] Keller | Generalized bounded variation and applications to piecewise monotonic transformations[END_REF] for one-dimensional dynamics, we refer the reader to [START_REF] Saussol | Absolutely continuous invariant measures for multidimensional expanding maps[END_REF] for a detailed presentation of that space in higher dimensions, as well as for the proof of its main properties. In particular, using the same notation as in [START_REF] Saussol | Absolutely continuous invariant measures for multidimensional expanding maps[END_REF], we use the following notion of variation:

var(f ) = sup 0<ε≤ε 0 ε -α R n osc(f, B ε (x)) dx, where osc(f, B ε (x)) = esssup Bε(x) f -essinf Bε(x) f.
In [START_REF] Saussol | Absolutely continuous invariant measures for multidimensional expanding maps[END_REF] it is proved that this notion of variation satisfies (V1)-( V3) and (V5)-(V7) and noted that (V4) is a consequence of a result in [START_REF] Keller | Generalized bounded variation and applications to piecewise monotonic transformations[END_REF]. We prove here that (V8) also holds. Assume that essinf f ≥ c and observe that

osc(1/f, B ε (x)) = esssup Bε(x) (1/f ) -essinf Bε(x) (1/f ) = 1/ essinf Bε(x) (f ) -1/ esssup Bε(x) (f ) = esssup Bε(x) (f ) -essinf Bε(x) (f ) (esssup Bε(x) (f ))(essinf Bε(x) (f )) ≤ osc(f, B ε (x)) (essinf Bε(x) (f )) 2 ≤ 1 c 2 osc(f, B ε (x)), which readily implies that var(1/f ) ≤ var(f ) (essinf f ) 2 .
One can consider piecewise C 2 expanding maps on compact subdomains of R d with C 2 boundary. As in the one-dimensional situation, condition (H3) holds because N is finite. In the local setting the conditions (H4) and (H5) can be verified as in the one-dimensional case using the results of [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF] (see also Th. 7.7 in [START_REF] Haydn | Almost sure invariance principle for sequential and non-stationary dynamical systems[END_REF]). In the nonlocal setting, for maps of the type considered in [3, Appendix B], it is likely that (H4) and (H5) can be obtained as in the one-dimensional case.

2.3.

Further properties of the random ACIM. Let µ ω be the measure on X given by dµ ω = h ω dm for ω ∈ Ω. We have the following important consequence of (H4), which establishes the appropriate decay of correlations result that will be used later on.

Lemma 1. There exists K > 0 and ρ ∈ (0, 1) such that

L n ω (φh ω )ψ dm -φ dµ ω • ψ dµ σ n (ω) ≤ Kρ n ψ ∞ • φ var , (5)
for n ≥ 0, ψ ∈ L ∞ (X, m) and φ ∈ BV (X, m).

Proof. We consider two cases. Assume first that φ dµ ω = φh ω dm = 0. Then, it follows from (H4) that

L n ω (φh ω )ψ dm -φ dµ ω • ψ dµ σ n ω = L n ω (φh ω )ψ dm ≤ ψ ∞ • L n ω (φh ω ) 1 ≤ ψ ∞ • L n ω (φh ω ) BV ≤ Ke -λn φh ω BV • ψ ∞ , and
thus (5) follows from (1) and ( 3). Now we consider the case when φ dµ ω = 0. We have

L n ω (φh ω )ψ dm -φ dµ ω • ψ dµ σ n (ω) = L n ω (φh ω )ψ dm -φh ω dm • ψh σ n (ω) dm ≤ ψ ∞ • L n ω (φh ω ) - φh ω dm h σ n (ω) dm = ψ ∞ • φh ω dm • L n ω (Φ -h ω ) dm ≤ ψ ∞ • φh ω dm • L n ω (Φ -h ω ) BV ,
where

φh ω = φh ω dm Φ.
Note that (Φ -h ω ) dm = 0 and thus using (H4),

ψ ∞ • φh ω dm • L n ω (Φ -h ω ) BV ≤ Ke -λn ψ ∞ • φh ω dm • Φ -h ω BV ≤ Ke -λn ψ ∞ • φ -φh ω dm h ω BV .
Hence, it follows from (1) and (3) that

L n ω (φh ω )ψ dm -φ dµ ω • ψ dµ σ n (ω) ≤ K e -λn ψ ∞ • φ BV
for some K > 0 and thus (5) follows from the observation that • BV ≤ • var .

Remark 1. We would like to emphasize that (5) is a special case of a more general decay of correlation result obtained in [START_REF] Buzzi | Exponential Decay of Correlations for Random Lasota-Yorke Maps[END_REF] which does not require (H4) and yields (5) but with

K = K(ω).
Finally, we prove that condition (H5) implies that we have a uniform lower bound for h ω .

Lemma 2. We have that essinf h ω ≥ c/2, for a.e. ω ∈ Ω.

Proof. We note that

h ω = L n σ -n (ω) 1 X -(L n σ -n (ω) 1 X -h ω ) and thus it follows from (H5) that essinf h ω ≥ c -L n σ -n (ω) 1 -h ω ) ∞ ≥ c -C var L n σ -n (ω) 1 -h ω BV . (7)
On the other hand, it follows from (H4) and (3) that

L n σ -n (ω) 1 -h ω BV = L n σ -n (ω) (1 -h σ -n (ω) ) BV ≤ Ke -λn 1 -h σ -n (ω) BV ≤ Ke -λn ,
for some K > 0. Choosing n such that C var Ke -λn ≤ c/2, it follows from ( 7) that ( 6) holds.

Remark 2. We now briefly compare our setting with that in [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF]. In the latter, the space X is replaced by a foliation Ξ ω := {ξ ∈ Ξ : (ξ, ω) ∈ Ξ}, where Ξ is some measurable space and ω belongs to the base space Ω. On the fibered subset Ξ it is defined the skew map τ (ξ, ω) = (f ω ξ, σω), with the associated fiber maps f ω : Ξ ω → Ξ σω . In our situation the Ξ ω 's for all ω coincide with the set X and all f ω : X → X are endomorphisms of X with some regularity property; moreover the previous skew transformation will still hold on the product space X × Ω, see [START_REF] Arnold | Random Dynamical Systems[END_REF]. Consequently, the conformal measure m is also allowed to depend on ω. In principle, all the arguments in the present paper could also be extended to this more general setting. However, we refrain from doing so since it would require heavy notation and more importantly since it is hard to verify that conditions like (H4) and (H5) hold in this more general setting unless fibers Ξ ω can be identified in some natural way (like for example tangent spaces at different points on a manifold).

2.4.

Statement of the main results. We are now ready to state our main result. We will consider an observable ψ : Ω × X → R. Let ψ ω = ψ(ω, •), ω ∈ Ω and assume that

sup ω∈Ω ψ ω BV < ∞. (8)
We also introduce centered observable

ψω = ψ ω -ψ ω dµ ω , ω ∈ Ω.
and we consider the associated Birkhoff sum

n-1 k=0 ψσ k ω • f k ω , and the variance τ 2 n = E ω n-1 k=0 ψσ k ω • f k ω 2
. The Almost Sure Invariance Principle is a matching of the trajectories of the dynamical system with a Brownian motion in such a way that the error is negligible in comparison with the Birkhoff sum. Limit theorems such as the central limit theorem, the functional central limit theorem and the law of the iterated logarithm transfer from the Brownian motion to time-series generated by observations on the dynamical system: these last results will therefore be immediate consequences of our proof of the ASIP for random Lasota-Yorke maps : Theorem 1. Let us consider the family of uniformly good random Lasota-Yorke maps. Then the variance τ 2 n will grow linearly as τ 2 n ∼ nΣ 2 and two cases will present: (i) either Σ = 0, and this is equivalent to the existence of φ ∈ L 2 (Ω × X) such that (co-boundary condition)

ψ = φ -φ • τ. (9)
(ii) or Σ 2 > 0 and in this case for P-a.e. ω ∈ Ω and d ∈ (0, 1/2), by enlarging probability space (X, B, µ ω ) if necessary, it is possible to find a sequence (Z k ) k of independent centered Gaussian random variables such that

sup 1≤k≤n n k=1 ( ψσ k ω • f k ω ) - n k=1 Z k = o((n 1/2+d (|log n 1/2-d | + log log n 1/2+d
)) 1/2 ), a.s.

Reverse martingale construction

In this section we construct the reverse martingale (or the reverse martingale difference) and establish various useful estimates that will play an important role in the rest of the paper. Indeed, the proof of our main result (Theorem 1) will be obtained as a consequence of the recent result by Cuny and Merlevède [START_REF] Cuny | Strong invariance principles with rate for "reverse" martingales and applications[END_REF] applied to our reverse martingale.

For ω ∈ Ω and k ∈ N, let T k ω := (f k ω ) -1 (B). Furthermore, for a measurable map φ : X → R and a σ-algebra H on X, E ω (φ|H) will denote the conditional expectation of φ with respect to H and the measure µ ω . Moreover, from now on, for k ∈ N and ω ∈ Ω we will write σ k ω instead of σ k (ω). We begin by the following technical lemma.

Lemma 3. We have

E ω (φ • f l ω |T n ω ) = L n-l σ l ω (h σ l ω φ) h σ n ω • f n ω , (10)
for each ω ∈ Ω and 0 ≤ l ≤ n.

Proof. We note that the right-hand side of ( 10) is measurable with respect to T n ω . Take now an arbitrary A ∈ T n ω and write it in the form A = (f n ω ) -1 (B) for some B ∈ B. We have

A φ • f l ω dµ ω = X (φ • f l ω )1 A dµ ω = X (φ • f l ω ) • (1 B • f n ω ) dµ ω = X φ(1 B • f n-l σ l ω ) dµ σ l ω = X h σ l ω φ(1 B • f n-l σ l ω ) dm = X L n-l σ l ω (h σ l ω φ)1 B dm = X L n-l σ l ω (h σ l ω φ) h σ n ω 1 B dµ σ n ω = X L n-l σ l ω (h σ l ω φ) h σ n ω • f n ω (1 B • f n ω ) dµ ω = X L n-l σ l ω (h σ l ω φ) h σ n ω • f n ω 1 A dµ ω = A L n-l σ l ω (h σ l ω φ) h σ n ω • f n ω dµ ω ,
which proves [START_REF] Keller | Generalized bounded variation and applications to piecewise monotonic transformations[END_REF].

We now return to the observable ψ ω introduced in (8) and its centered companion

ψω = ψ ω -ψ ω dµ ω , ω ∈ Ω. Set M n = ψσ n ω + G n -G n+1 • f σ n ω , n ≥ 0, (11)
where G 0 = 0 and

G k+1 = L σ k ω ( ψσ k ω h σ k ω + G k h σ k ω ) h σ k+1 ω , k ≥ 0. ( 12
)
We emphasize that M n and G n depend on ω. However, in order to avoid complicating the notation, we will not make this dependence explicit. In preparation for the next proposition we need the following elementary result.

Lemma 4. We have

L ω ((ψ • f ω )φ) = ψL ω φ, for φ ∈ L 1 (X, m) and ψ ∈ L ∞ (X, m).
Proof. We first note that (ψ • f ω )φ ∈ L 1 (X, m). Moreover, for an arbitrary g ∈ L ∞ (X, m) we have that

X L ω ((ψ • f ω )φ)g dm = X (ψ • f ω )φ(g • f ω ) dm = X ((ψg) • f ω )φ dm = X (ψL ω φ)g dm,
which immediately implies the conclusion of the lemma. Now we prove that the sequence (M n • f n ω ) n is a reversed martingale (or the reversed martingale difference) with respect to the sequence of σ-algebras (T n ω ) n . Proposition 2. We have

E ω (M n • f n ω |T n+1 ω ) = 0.
Proof. It follows from Lemma 3 that

E ω (M n • f n ω |T n+1 ω ) = L σ n ω (h σ n ω M n ) h σ n+1 ω • f n+1 ω . ( 13
)
Moreover, by [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF] we have

L σ n ω (h σ n ω M n ) h σ n+1 ω = L σ n ω (h σ n ω ψσ n ω + h σ n ω G n -h σ n ω (G n+1 • f σ n ω )) h σ n+1 ω .
By Lemma 4,

L σ n ω (h σ n ω (G n+1 • f σ n ω )) = G n+1 L σ n ω h σ n ω = G n+1 h σ n+1 ω ,
and thus it follows from ( 12) that

L σ n ω (h σ n ω M n ) h σ n+1 ω = 0.
This conclusion of the lemma now follows readily from [START_REF] Philip | Almost sure invariance principle for sums of weakly dependent random variables[END_REF].

We now establish several auxiliary results that will be used in the following section. These results estimate various norms of functions related to M n and G n , defined in [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF] and [START_REF] Liverani | Decay of correlations for piecewise expanding maps[END_REF], respectively. Lemma 5. We have that sup

n≥0 G n BV < ∞.
Proof. By iterating [START_REF] Liverani | Decay of correlations for piecewise expanding maps[END_REF], we obtain

G n = 1 h σ n ω n-1 j=0 L n-j σ j ω ( ψσ j ω h σ j ω ), n ∈ N. (14)
We note that ψσ j ω h σ j ω dm = ψσ j ω dµ σ j ω = 0,

and thus it follows from (H4) that

n-1 j=0 L n-j σ j ω ( ψσ j ω h σ j ω ) BV ≤ K n-1 j=0
e -λ(n-j) ψσ j ω h σ j ω BV , for each n ∈ N which together with (V8), ( 1), ( 3), ( 6) and ( 8) implies the conclusion of the lemma.

Lemma 6. We have that

sup n≥0 M 2 n BV < ∞.
Proof. In view of ( 8), [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF] and Lemma 5, it is sufficient to show that

sup n≥0 G n+1 • f σ n ω BV < ∞.
However, this follows directly from (H3) and Lemma 5.

Lemma 7. We have that

sup n≥0 E ω (M 2 n • f n ω |T n+1 ω ) ∞ < ∞.
Proof. It follows from Lemma 3 that

E ω (M 2 n • f n ω |T n+1 ω ) = L σ n ω (h σ n ω M 2 n ) h σ n+1 ω • f n+1 ω ,
and thus, recalling [START_REF] Cuny | Strong invariance principles with rate for "reverse" martingales and applications[END_REF],

sup n≥0 E ω (M 2 n • f n ω |T n+1 ω ) ∞ ≤ 1 c L σ n ω (h σ n ω M 2 n ) ∞ .
Taking into account ( 1), (H2), (3), Lemma 6 and the fact that

• ∞ ≤ C var • BV (see ( V3
)) we obtain the conclusion of the lemma.

Sprindzuk's Theorem and consequences

The main tool in establishing the almost sure invariance principle is the recent result by Cuny and Merlevède (quoted in our Theorem 3 in Section 5). However, in order to verify the assumptions of that theorem we will first need to apply the following classical result due to Sprindzuk [START_REF] Sprindzuk | Metric theory of Diophantine approximations[END_REF]. Theorem 2. Let (Ω, B, µ) be a probability space and let (f k ) k be a sequence of nonnegative and measurable functions on Ω. Moreover, let (g k ) k and (h k ) k be bounded sequences of real numbers such that 0 ≤ g k ≤ h k . Assume that there exists C > 0 such that

m<k≤n (f k (ω) -g k ) 2 dµ ≤ C m<k≤n h k (16)
for m, n ∈ N such that m < n. Then, for every ε > 0

1≤k≤n f k (ω) = 1≤k≤n g k + O(Θ 1/2 (n) log 3/2+ε Θ(n)),
for µ-a.e. ω ∈ Ω, where

Θ(n) = 1≤k≤n h k . Lemma 8. For each ε > 0, n-1 k=0 E ω (M 2 k • f k ω |T k+1 ω ) = n-1 k=0 E ω (M 2 k • f k ω ) + O(Θ 1/2 (n) log 3/2+ε Θ(n)),
for µ-a.e. ω ∈ Ω, where

Θ(n) = n-1 k=0 (E ω (M 2 k • f k ω ) + M 2 k var ). (17)
Proof. Fix ω ∈ Ω. We want to apply Theorem 2 to

f k = E ω (M 2 k • f k ω |T k+1 ω ) and g k = E ω (M 2 k • f k ω ).
We have that m<k≤n

E ω (M 2 k • f k ω |T k+1 ω ) - m<k≤n E ω (M 2 k • f k ω ) 2 dµ ω = m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) 2 dµ ω -2 m<k≤n E ω (M 2 k • f k ω ) m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) dµ ω + m<k≤n E ω (M 2 k • f k ω ) 2 = m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) 2 dµ ω - m<k≤n E ω (M 2 k • f k ω ) 2 = m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) 2 dµ ω + 2 m<i<j≤n E ω (M 2 i • f i ω |T i+1 ω ) • E ω (M 2 j • f j ω |T j+1 ω ) dµ ω - m<k≤n E ω (M 2 k • f k ω ) 2 -2 m<i<j≤n E ω (M 2 i • f i ω ) • E ω (M 2 j • f j ω ) ≤ m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) 2 dµ ω + 2 m<i<j≤n E ω (M 2 i • f i ω |T i+1 ω ) • E ω (M 2 j • f j ω |T j+1 ω ) dµ ω -2 m<i<j≤n E ω (M 2 i • f i ω ) • E ω (M 2 j • f j ω ). (18)
On the other hand, it follows from Lemma 3 that for i < j we have

E ω (M 2 i • f i ω |T i+1 ω ) • E ω (M 2 j • f j ω |T j+1 ω ) dµ ω = L σ i ω (h σ i ω M 2 i ) h σ i+1 ω • f i+1 ω • L σ j ω (h σ j ω M 2 j ) h σ j+1 ω • f j+1 ω dµ ω = L σ i ω (h σ i ω M 2 i ) h σ i+1 ω • L σ j ω (h σ j ω M 2 j ) h σ j+1 ω • f j-i σ i+1 ω dµ σ i+1 ω = L σ i ω (h σ i ω M 2 i ) • L σ j ω (h σ j ω M 2 j ) h σ j+1 ω • f j-i σ i+1 ω dm = L j-i+1 σ i ω (h σ i ω M 2 i ) • L σ j ω (h σ j ω M 2 j ) h σ j+1 ω dm.
Moreover,

E ω (M 2 i • f i ω ) = M 2 i dµ σ i ω and E ω (M 2 j • f j ω ) = (M 2 j • f j ω ) dµ ω = M 2 j dµ σ j ω = M 2 j h σ j ω dm = L σ j ω (M 2 j h σ j ω ) dm = L σ j ω (M 2 j h σ j ω ) h σ j+1 ω dµ σ j+1 ω .
Hence,

E ω (M 2 i • f i ω |T i+1 ω ) • E ω (M 2 j • f j ω |T j+1 ω ) dµ ω -E ω (M 2 i • f i ω ) • E ω (M 2 j • f j ω ) = L j-i+1 σ i ω (h σ i ω M 2 i ) • L σ j ω (h σ j ω M 2 j ) h σ j+1 ω dm -M 2 i dµ σ i ω • L σ j ω (M 2 j h σ j ω ) h σ j+1 ω dµ σ j+1 ω .
Therefore, it follows from Lemma 1 that

E ω (M 2 i • f i ω |T i+1 ω ) • E ω (M 2 j • f j ω |T j+1 ω ) dµ ω -E ω (M 2 i • f i ω ) • E ω (M 2 j • f j ω ) ≤ Kρ j-i+1 L σ j ω (M 2 j h σ j ω ) h σ j+1 ω ∞ • M 2 i var .
Furthermore,

E ω (M 2 k • f k ω |T k+1 ω ) 2 dµ ω ≤ E ω (M 2 k • f k ω |T k+1 ω ) ∞ • E ω (M 2 k • f k ω ).
Thus, the last two inequalities combined with (18) imply that m<k≤n

E ω (M 2 k • f k ω |T k+1 ω ) - m<k≤n E ω (M 2 k • f k ω ) 2 dµ ω ≤ m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) 2 dµ ω + 2 m<i<j≤n E ω (M 2 i • f i ω |T i+1 ω ) • E ω (M 2 j • f j ω |T j+1 ω ) dµ ω -2 m<i<j≤n E ω (M 2 i • f i ω ) • E ω (M 2 j • f j ω ). ≤ m<k≤n E ω (M 2 k • f k ω |T k+1 ω ) ∞ • E ω (M 2 k • f k ω ) + 2K m<i<j≤n ρ j-i+1 L σ j ω (M 2 j h σ j ω ) h σ j+1 ω ∞ • M 2 i var ,
which combined with (H2), ( 3), [START_REF] Cuny | Strong invariance principles with rate for "reverse" martingales and applications[END_REF] and Lemmas 6 and 7 implies that (16) holds with

h k = E ω (M 2 k • f k ω ) + M 2 k var .
The conclusion of the lemma now follows directly from Theorem 2.

Proof of Theorem 1

The goal of this section is to establish the almost sure invariance principle by proving Theorem 1. It is based on the following theorem due to Cuny and Merlevède. Theorem 3 ([6]). Let (X n ) n be a sequence of square integrable random variables adapted to a non-increasing filtration (G n ) n . Assume that E(X n |G n+1 ) = 0 a.s., σ 2 n := n k=1 E(X 2 k ) → ∞ when n → ∞ and that sup n E(X 2 n ) < ∞. Moreover, let (a n ) n be a non-decreasing sequence of positive numbers such that the sequence (a n /σ 2 n ) n is non-increasing, (a n /σ n ) is non-decreasing and such that :

1. n k=1 (E(X 2 k |G k+1 ) -E(X 2 k )) = o(a n ) a.s.; (19) 2.
n≥1

a -v n E(|X n | 2v ) < ∞ for some 1 ≤ v ≤ 2. ( 20
)
Then, enlarging our probability space if necessary, it is possible to find a sequence

(Z k ) k of independent centered Gaussian variables with E(X 2 k ) = E(Z 2 k ) such that sup 1≤k≤n n k=1 X k - n k=1 Z k = o((a n (|log(σ 2 n /a n )| + log log a n )) 1/2 ), a.s.
In order to obtain the almost sure invariance principle for the sequence ( ψθ k ω • f k ω ) k , k ∈ N we will first apply Theorem 3 for

X n = M n • f n ω and G n = T n ω . We note that it follows from Lemma 8 that n k=1 (E(X 2 k |G k+1 ) -E(X 2 k )) = O(b n ), with b n = Θ 1/2 (n) log 3/2+ε Θ(n), (21)
and where Θ(n) is given by (17). On the other hand, it follows from Lemma 6 that Θ(n) ≤ Dn for some D > 0 and every n ∈ N and therefore (19) holds with

a n = n 1/2+d , (22)
for any d > 0. From now on, we take d ∈ (0, 1/2).

Lemma 9. There exists Σ 2 ≥ 0 such that

lim n→∞ 1 n E ω n-1 k=0 ψσ k ω • f k ω 2 = Σ 2 , for a.e. ω ∈ Ω. (23)
Proof. Note that

E ω n-1 k=0 ψσ k ω • f k ω 2 = n-1 k=0 E ω ( ψ2 σ k ω • f k ω ) + 2 0≤i<j≤n-1 E ω (( ψσ i ω • f i ω )( ψσ j ω • f j ω )) = n-1 k=0 E ω ( ψ2 σ k ω • f k ω ) + 2 n-1 i=0 n-1 j=i+1 E σ i ω ( ψσ i ω ( ψσ j ω • f j-i σ i ω )).
Using the skew product transformation τ from (2), it follows from Birkhoff's ergodic theorem that

lim n→∞ 1 n n-1 k=0 E ω ( ψ2 σ k ω • f k ω ) = lim n→∞ 1 n n-1 k=0 X ψ(τ k (ω, x)) 2 dµ ω (x) = Ω X ψ(ω, x) 2 dµ ω (x) dP(ω) = Ω×X ψ(ω, x) 2 dµ(ω, x),
for a.e. ω ∈ Ω, where µ is an invariant measure for τ given by

µ(A × B) = A µ ω (B) dP(ω), for A ∈ F and B ∈ B.
Furthermore, set

Ψ(ω) = ∞ n=1 X ψ(ω, x) ψ(τ n (ω, x)) dµ ω (x) = ∞ n=1 X L n ω ( ψω h ω ) ψσ n ω dm.
By ( 5) and ( 8), we have

|Ψ(ω)| ≤ ∞ n=1 X L n ω ( ψω h ω ) ψσ n ω dm ≤ K ∞ n=1 ρ n = Kρ 1 -ρ ,
for some K > 0 and a.e. ω ∈ Ω. In particular, Ψ ∈ L 1 (Ω) and thus it follows again from Birkhoff's ergodic theorem that

lim n→∞ 1 n n-1 i=0 Ψ(σ i ω) = Ω Ψ(ω) dP(ω) = ∞ n=1 Ω×X ψ(ω, x) ψ(τ n (ω, x)) dµ(ω, x), (24)
for a.e. ω ∈ Ω. In order to complete the proof of the lemma, we are going to show that

lim n→∞ 1 n n-1 i=0 n-1 j=i+1 E σ i ω ( ψσ i ω ( ψσ j ω • f j-i σ i ω )) - n-1 i=0 Ψ(σ i ω) = 0, (25)
for a.e. ω ∈ Ω. Using (5), we have that for a.e. ω ∈ Ω,

n-1 i=0 n-1 j=i+1 E σ i ω ( ψσ i ω ( ψσ j ω • f j-i σ i ω )) - n-1 i=0 Ψ(σ i ω) = n-1 i=0 n-1 j=i+1 E σ i ω ( ψσ i ω ( ψσ j ω • f j-i σ i ω )) - n-1 i=0 ∞ k=1 E σ i ω ( ψσ i ω ( ψσ k+i ω • f k σ i ω )) ≤ n-1 i=0 ∞ k=n-i E σ i ω ( ψσ i ω ( ψσ k+i ω • f k σ i ω )) = n-1 i=0 ∞ k=n-i X L k σ i ω ( ψσ i ω h σ i ω ) ψσ k+i ω dm ≤ K n-1 i=0 ∞ k=n-i ρ k = K ρ (1 -ρ) 2 ,
which readily implies (25). It follows from ( 24) and (25) that

lim n→∞ 1 n n-1 i=0 n-1 j=i+1 E σ i ω ( ψσ i ω ( ψσ j ω • f j-i σ i ω )) = ∞ n=1 Ω×X ψ(ω, x) ψ(τ n (ω, x)) dµ(ω, x)
for a.e. ω ∈ Ω and therefore (23) holds with

Σ 2 = Ω×X ψ(ω, x) 2 dµ(ω, x) + 2 ∞ n=1 Ω×X ψ(ω, x) ψ(τ n (ω, x)) dµ(ω, x). (26)
Finally, we note that it follows readily from (23) that Σ 2 ≥ 0 and the proof of the lemma is completed.

We now present necessary and sufficient conditions under which Σ 2 = 0. We note that a similar result is stated in [11, (2.10)] with ψ • τ instead of ψ in (27).

Proposition 3. We have that Σ 2 = 0 if and only if there exists φ ∈ L 2 (Ω × X) such that

ψ = φ -φ • τ. ( 27
)
Proof. We first observe that Assume now that Σ 2 = 0. Then, it follows from the above equality and (26) that On the other hand, by [START_REF] Conze | Limit theorems for sequential expanding dynamical systems on [0, 1], In: Ergodic Theory and Related Fields[END_REF] we have that Ω×X ψ( ψ • τ k ) dµ → 0 exponentially fast when k → ∞ and hence, it follows from (28) that the sequence (X n ) n defined by

X n (ω, x) = n-1 k=0 ψ(τ k (ω, x)), ω ∈ Ω, x ∈ X
is bounded in L 2 (Ω × X). Thus, it has a subsequence (X n k ) k which converges weakly to some φ ∈ L 2 (Ω × X). We claim that φ satisfies (27). Indeed, take an arbitrary g = 1 A×B , where A ∈ F and B ∈ B and observe that g ∈ L 2 (Ω × X) and

Ω×X g( ψ -φ + φ • τ ) = lim k→∞ Ω×X g( ψ -X n k + X n k • τ ) dµ = lim k→∞ Ω×X g( ψ • τ n k ) dµ = 0,
where in the last equality we used (5) again. Therefore, ψ -φ + φ • τ = 0 which readily implies (27). Suppose now that there exists φ ∈ L 2 (Ω × X) satisfying (27). Then,

1 √ n n-1 k=0 ψ • τ k = 1 √ n (φ -φ • τ n ),
and thus

1 √ n n-1 k=0 ψ • τ k L 2 (Ω×X) ≤ 2 √ n φ L 2 (Ω×X) → 0,
when n → ∞. Therefore, it follows by integrating (23) over Ω that

Σ 2 = lim n→∞ 1 √ n n-1 k=0 ψ • τ k 2 L 2 (Ω×X) = 0.
This concludes the proof of the proposition.

In the rest of the paper we assume that Σ 2 > 0. We also need the following lemmas.

Lemma 10. We have that E ω (X i X j ) = 0, for i < j.

Proof. By Lemma 2, we conclude that E ω (M i • f i ω |T i+1 ω ) = 0. Moreover, we note that M j • f j ω is measurable with respect to T i+1 ω and thus

E ω ((M j • f j ω )(M i • f i ω )|T i+1 ω ) = (M j • f j ω )E ω (M i • f i ω |T i+1 ω ) = 0.
The conclusion of the lemma now follows simply by integrating the above equality.

In what follows, we write a n ∼ b n if there exists c ∈ R\{0} such that lim n→∞ a n /b n = c.

Lemma 11. We have that σ 2 n → ∞ as n → ∞. Proof. It follows from [START_REF] Kifer | Limit theorems for random transformations and processes in random environments[END_REF] that

n-1 k=0 X k = n-1 k=0 ψσ k ω • f k ω -G n • f n ω , (29)
and thus,

n-1 k=0 X k 2 = n-1 k=0 ψσ k ω • f k ω 2 -2(G n • f n ω ) n-1 k=0 ψσ k ω • f k ω + (G 2 n • f n ω ). (30)
By Lemma 9 and the assumption Σ 2 > 0,

τ 2 n := E ω n-1 k=0 ψσ k ω • f k ω 2 → ∞. (31)

  k (ω, x)) dµ(ω, x) + 2 k (ω, x) ψ(τ j (ω, x)) dµ(ω, x) x) ψ(τ k-j (ω, x)) dµ(ω, x) x) ψ(τ k (ω, x)) dµ(ω, x), , x) ψ(τ k (ω, x)) dµ(ω, x) , x) ψ(τ k (ω, x)) dµ(ω, x).

Ω×X n- 1 k=0ψ • τ k 2 dµ

 12 • τ k ) dµ. (28)

Acknowledgments

The research of DD and GF was supported by the Australian Research Council Discovery Project DP150100017 and, in part, by the Croatian Science Foundation under the project IP-2014-09-228. CGT acknowledges the support of Australian Research Council Discovery Early Career Researcher Award (DECRA) DE160100147. SV was supported by the ANR-Project Perturbations and by the Leverhulme Trust for support thorough the Network Grant IN-2014-021. SV and CGT warmly thank the School of Mathematics and Statistics of the University of New South Wales for the kind hospitality during the preparation of this work. SV acknowledges several discussions with M. Abdelkader and

On the other hand, it follows from ( 8), (30) and Lemma 5 that

By Lemma 10 and (32), we have that

which together with (31) implies the desired conclusion of Lemma 11.

Lemma 12. There exists n 0 ∈ N such that the sequence (a n /σ 2 n ) n≥n 0 is non-increasing and that the sequence (a n /σ n ) n≥n 0 is non-decreasing.

Proof. It follows from Lemma 9 and (33) that

Therefore (22) yields,

and the conclusion of the lemma follows readily from the assumption that d ∈ (0, 1/2).

Since the conclusion of Theorem 3 concerns the tails of (a n ) n and (σ n ) n , it will remain valid if the monotonicity hypotheses for (a n /σ 2 n ) n and (a n /σ n ) n hold for sufficiently large n, and those are verified in Lemma 12. Finally, we show that (20) holds with v = 2.

Lemma 13. We have that

Now we can conclude the proof or our main result.

Proof of Theorem 1. Using Theorem 3, we obtain the almost sure invariance principle for the sequence

The almost sure invariance principle for the sequence


Yoann Thomas

Fred Jean

Laure Pecquerie

Gonçalo Marques

Sébastien Hervé

Valérie Dantec

Jonathan Flye-Sainte-Marie

DEBSe@: an e-learning platform to introduce the Dynamic Energy Budget (DEB) theory for marine ecology, fisheries sciences and aquaculture applications

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Prac%cal informa%on

Course schedule: two sessions of 4 to 5 weeks will be open every years for the online course, followed by a class :me of two weeks. Target audience: Master and PhD students & researchers who might be interested in an introductory course on DEB theory. Contacts: For more informa:on, please contact jonathan.flye@univ--brest.fr ; laure.pecquerie@ird.fr

Content of the online course

Learning modules Serious--game "Béber project" A virtual field experience

Concept:

Rearing individuals from a cohort of DEBSe@ organisms (see Figures) that will grow, develop and reproduce following the rules of the standard DEB model, according to the environmental condi:ons defined by students. Each student will create its own data set and then will use it:

• for the first applica:on of the DEB model build with R.

• to learn how to es:mate the parameters using the add--my--pet procedure a and the DEBtool rou:nes, recently translated in R.

Objec%ves:

• To s:mulate the course interac:vity.

• To create a personal dataset (environmental and biometric :me--series).

Context

DEBSe@ provides an introductory course on Dynamic Energy Budget (DEB) theory to master and PhD students in marine biological sciences worldwide.

DEBSe@ follows a flipped classroom concept, with:

• an interna%onal Small Private Online Course (SPOC) followed by

• a class in presence of a teacher to apply what have been learned and to go to higher--order tasks

The online course will be available in French, Spanish and English from the MOODLE plaMorm managed by the Université de Bretagne Occidentale (UBO, Brest, France).

DEBSe@ is complementary to the Massive Open Online DEB--Course (DEB--MOOC a , Akvaplan--Niva, Norway and IST, Portugal), targe:ng advanced researchers and PhD students.

Life cycle of the DEBSe@ "Béber" organism, following DEB theory concepts

: an e--learning plaLorm to introduce the Dynamic Energy Budget (DEB) theory for marine ecology, fisheries sciences and aquaculture applica:ons Yoann Thomas 1 , Laure Pecquerie 2 , Frédéric Jean 1 , Gonçalo Marques 3 , Sébas:en Hervé 1 , Valérie Dantec 4 , Jonathan Flye--Sainte--Marie 1 3D--animaBons will be available to illustrate the various state and processes (feeding, growth etc.) of the DEBSe@ "Béber" organism

-

  -Online course plaLorm: MOODLE --Access area tailored to individual student --Auto--assessment --Several asset--based materials (graphics, anima:ons etc.) --Interac:ve ac:vi:es (see below) --Service--based features (online assessments, forums, news etc.) --External links (Web pages, books, published ar:cles etc.) --Feedback on the course enabling its improvement

1

  Université de Bretagne Occidentale, UMR LEMAR --6539 CNRS/UBO/IRD/Ifremer, France 2 Ins:tut de Recherche pour le Développement, UMR LEMAR --6539 CNRS/UBO/IRD/Ifremer, France 3 Ins:tuto Superior Técnico, Universidade Técnica de Lisboa, Portugal 4 Université de Bretagne Occidentale, SIAME, France a hip://www.debtheory.org/ Conceptual scheme of the DEB model under construcBon, illustraBng one of the learning module of the course the DEB model and a theore:cal data set generated during the serious--game "Béber project"


Eugene Belilovsky
email: eugene.belilovsky@inria.fr

Gael Varoquaux
email: gael.varoquaux@inria.fr

Matthew Blaschko
email: matthew.blaschko@esat.kuleuven.be

Testing for Differences in Gaussian Graphical Models: Applications to Brain Connectivity

Functional brain networks are well described and estimated from data with Gaussian Graphical Models (GGMs), e.g. using sparse inverse covariance estimators. Comparing functional connectivity of subjects in two populations calls for comparing these estimated GGMs. Our goal is to identify differences in GGMs known to have similar structure. We characterize the uncertainty of differences with confidence intervals obtained using a parametric distribution on parameters of a sparse estimator. Sparse penalties enable statistical guarantees and interpretable models even in high-dimensional and low-sample settings. Characterizing the distributions of sparse models is inherently challenging as the penalties produce a biased estimator. Recent work invokes the sparsity assumptions to effectively remove the bias from a sparse estimator such as the lasso. These distributions can be used to give confidence intervals on edges in GGMs, and by extension their differences. However, in the case of comparing GGMs, these estimators do not make use of any assumed joint structure among the GGMs. Inspired by priors from brain functional connectivity we derive the distribution of parameter differences under a joint penalty when parameters are known to be sparse in the difference. This leads us to introduce the debiased multi-task fused lasso, whose distribution can be characterized in an efficient manner. We then show how the debiased lasso and multi-task fused lasso can be used to obtain confidence intervals on edge differences in GGMs. We validate the techniques proposed on a set of synthetic examples as well as neuro-imaging dataset created for the study of autism.

Introduction

Gaussian graphical models describe well interactions in many real-world systems. For instance, correlations in brain activity reveal brain interactions between distant regions, a process know as functional connectivity. Functional connectivity is an interesting probe on brain mechanisms as it persists in the absence of tasks (the so-called "resting-state") and is thus applicable to study populations of impaired subjects, as in neurologic or psychiatric diseases [START_REF] Castellanos | Clinical applications of the functional connectome[END_REF]. From a formal standpoint, Gaussian graphical models are well suited to estimate brain connections from functional Magnetic Resonance Imaging (fMRI) signals [START_REF] Smith | Network modelling methods for fMRI[END_REF][START_REF] Varoquaux | Brain covariance selection: Better individual functional connectivity models using population prior[END_REF]. A set of brain regions and related functional connections is then called a functional connectome [START_REF] Varoquaux | Learning and comparing functional connectomes across subjects[END_REF][START_REF] Castellanos | Clinical applications of the functional connectome[END_REF]. Its variation across subjects can capture cognition [START_REF] Richiardi | Decoding brain states from fMRI connectivity graphs[END_REF][START_REF] Shirer | Decoding subject-driven cognitive states with whole-brain connectivity patterns[END_REF] or pathology [START_REF] Kelly | Characterizing variation in the functional connectome: Promise and pitfalls[END_REF][START_REF] Castellanos | Clinical applications of the functional connectome[END_REF]. However, the effects of pathologies are often very small, as resting-state fMRI is a weakly-constrained and noisy imaging modality, and the number of subjects in a study is often small given the cost of imaging. Statistical power is then a major concern [START_REF] Button | Power failure: Why small sample size undermines the reliability of neuroscience[END_REF]. The statistical challenge is to increase the power to detect differences between Gaussian graphical models in the small-sample regime.

In these settings, estimation and comparison of Gaussian graphical models fall in the range of high-dimensional statistics: the number of degrees of freedom in the data is small compared to 30th Conference on Neural Information Processing Systems (NIPS 2016), Barcelona, Spain.

the dimensionality of the model. In this regime, sparsity-promoting 1 -based penalties can make estimation well-posed and recover good estimation performance despite the scarcity of data [START_REF] Tibshirani | Regression shrinkage and selection via the lasso[END_REF][START_REF] Friedman | Sparse inverse covariance estimation with the graphical lasso[END_REF][START_REF] Meinshausen | High-dimensional graphs and variable selection with the lasso[END_REF][START_REF] Danaher | The joint graphical lasso for inverse covariance estimation across multiple classes[END_REF][START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF]. These encompass sparse regression methods such as the lasso or recovery methods such as basis pursuit, and can be applied to estimation of Gaussian graphical models with approaches such as the graphical lasso [START_REF] Friedman | Sparse inverse covariance estimation with the graphical lasso[END_REF]. There is now a wide body of literature which demonstrates the statistical properties of these methods [START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF]. Crucial to applications in medicine or neuroscience, recent work characterizes the uncertainty, with confidence intervals and p-values, of the parameters selected by these methods [START_REF] Janková | Confidence intervals for high-dimensional inverse covariance estimation[END_REF][START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF][START_REF] Lockhart | A significance test for the lasso[END_REF][START_REF] Sell | Adaptive testing for the graphical lasso[END_REF]. These works focus primarily on the lasso and graphical lasso.

Approaches to estimate statistical significance on sparse models fall into several general categories: (a) non-parameteric sampling based methods which are inherently expensive and have difficult limiting distributions [START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF][START_REF] Narayan | Mixed effects models to find differences in multi-subject functional connectivity[END_REF][START_REF] Mota | Randomized parcellation based inference[END_REF], (b) characterizations of the distribution of new parameters that enter a model along a regularization path [START_REF] Lockhart | A significance test for the lasso[END_REF][START_REF] Sell | Adaptive testing for the graphical lasso[END_REF], or (c) for a particular regularization parameter, debiasing the solution to obtain a new consistent estimator with known distribution [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF][START_REF] Janková | Confidence intervals for high-dimensional inverse covariance estimation[END_REF][START_REF] Van De Geer | On asymptotically optimal confidence regions and tests for high-dimensional models[END_REF]. While some of the latter work has been used to characterize confidence intervals on network edge selection, there is no result, to our knowledge, on the important problem of identifying differences in networks. Here the confidence on the result is even more critical, as the differences are the direct outcome used for neuroscience research or medical practice, and it is important to provide the practitioner a measure of the uncertainty.

Here, we consider the setting of two datasets known to have very similar underlying signals, but which individually may not be very sparse. A motivating example is determining the difference in brain networks of subjects from different groups: population analysis of connectomes [START_REF] Varoquaux | Learning and comparing functional connectomes across subjects[END_REF][START_REF] Kelly | Characterizing variation in the functional connectome: Promise and pitfalls[END_REF]. Recent literature in neuroscience [START_REF] Markov | Cortical high-density counterstream architectures[END_REF] has suggested functional networks are not sparse. On the other hand, differences in connections across subjects should be sparse. Indeed the link between functional and anatomical brain networks [START_REF] Honey | Predicting human resting-state functional connectivity from structural connectivity[END_REF] suggests they should not differ drastically from one subject to another. From a neuroscientific standpoint we are interested in determining which edges between two populations (e.g. autistic and non-autistic) are different. Furthermore we want to provide confidence-intervals on our results. We particularly focus on the setting where one dataset is larger than the other. In many applications it is more difficult to collect one group (e.g. individuals with specific pathologies) than another.

We introduce an estimator tailored to this goal: the debiased multi-task fused lasso. We show that, when the underlying parameter differences are indeed sparse, we can obtain a tractable Gaussian distribution for the parameter difference. This closed-form distribution underpins accurate hypothesis testing and confidence intervals. We then use the relationship between nodewise regression and the inverse covariance matrix to apply our estimator to learning differences of Gaussian graphical models.

The paper is organized as follows. In Section 2 we review previous work on learning of GGMs and the debiased lasso. Section 3 discusses a joint debiasing procedure that specifically debiases the difference estimator. In Section 3.1 we introduce the debiased multi-task fused lasso and show how it can be used to learn parameter differences in linear models. In Section 3.2, we show how these results can be used for GGMs. In Section 4 we validate our approach on synthetic and fMRI data.

Background and Related Work

Debiased Lasso A central starting point for our work is the debiased lasso [START_REF] Van De Geer | On asymptotically optimal confidence regions and tests for high-dimensional models[END_REF][START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF]. Here one considers the linear regression model, Y = Xβ + , with data matrix X and output Y , corrupted by ∼ N (0, σ 2 I) noise. The lasso estimator is formulated as follows:

βλ = arg min β 1 n Y -Xβ 2 + λ β 1 (1)
The KKT conditions give kλ = 1 n X T (Y -Xβ), where k is the subgradient of λ β 1 . The debiased lasso estimator [START_REF] Van De Geer | On asymptotically optimal confidence regions and tests for high-dimensional models[END_REF][START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF] is then formulated as βλ u = βλ +M kλ for some M that is constructed to give guarantees on the asymptotic distribution of βλ u . Note that this estimator is not strictly unbiased in the finite sample case, but has a bias that rapidly approaches zero (w.r.t. n) if M is chosen appropriately, the true regressor β is indeed sparse, and the design matrix satistifes a certain restricted eigenvalue property [START_REF] Van De Geer | On asymptotically optimal confidence regions and tests for high-dimensional models[END_REF][START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF]. We decompose the difference of this debiased estimator and the truth as follows:

βλ u -β = 1 n M X T -(M Σ -I)( β -β) (2)
The first term is Gaussian and the second term is responsible for the bias. Using Holder's inequality the second term can be bounded by M Σ-I ∞ β -β 1 . The first part of which we can bound using an appropriate selection of M while the second part is bounded by our implicit sparsity assumptions coming from lasso theory [START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF]. Two approaches from the recent literature discuss how one can select M to appropriately debias this estimate. In [START_REF] Van De Geer | On asymptotically optimal confidence regions and tests for high-dimensional models[END_REF] it suffices to use nodewise regression to learn an inverse covariance matrix which guarantees constraints on M Σ -I ∞ . A second approach by [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF] proposes to solve a quadratic program to directly minimize the variance of the debiased estimator while constraining M Σ -I ∞ to induce sufficiently small bias.

Intuitively the construction of βλ u allows us to trade variance and bias via the M matrix. This allows us to overcome a naive bias-variance tradeoff by leveraging the sparsity assumptions that bound β -β 1 . In the sequel we expand this idea to the case of debiased parameter difference estimates and sparsity assumptions on the parameter differences.

In the context of GGMs, the debiased lasso can gives us an estimator that asymptotically converges to the partial correlations. As highlighted by [START_REF] Waldorp | Testing for graph differences using the desparsified lasso in high-dimensional data[END_REF] we can thus use the debiased lasso to obtain difference estimators with known distributions. This allows us to obtain confidence intervals on edge differences between Gaussian graphical models. We discuss this further in the sequel.

Gaussian Graphical Model Structure Learning A standard approach to estimating Gaussian graphical models in high dimensions is to assume sparsity of the precision matrix and have a constraint which limits the number of non-zero entries of the precision matrix. This constraint can be achieved with a 1 -norm regularizer as in the popular graphical lasso [START_REF] Friedman | Sparse inverse covariance estimation with the graphical lasso[END_REF]. Many variants of this approach that incorporate further structural assumptions have been proposed [START_REF] Honorio | Multi-task learning of Gaussian graphical models[END_REF][START_REF] Danaher | The joint graphical lasso for inverse covariance estimation across multiple classes[END_REF][START_REF] Mohan | Structured learning of Gaussian graphical models[END_REF].

An alternative solution to inducing sparsity on the precision matrix indirectly is neighborhood 1 regression from [START_REF] Meinshausen | High-dimensional graphs and variable selection with the lasso[END_REF]. Here the authors make use of a long known property that connects the entries of the precision matrix to the problem of regression of one variable on all the others [START_REF] Marsaglia | Conditional means and covariances of normal variables with singular covariance matrix[END_REF]. This property is critical to our proposed estimation as it allows relating regression models to finding edges connected to specific nodes in the GGM.

GGMs have been found to be good at recovering the main brain networks from fMRI data [START_REF] Smith | Network modelling methods for fMRI[END_REF][START_REF] Varoquaux | Brain covariance selection: Better individual functional connectivity models using population prior[END_REF]. Yet, recent work in neuroscience has showed that the structural wiring of the brain did not correspond to a very sparse network [START_REF] Markov | Cortical high-density counterstream architectures[END_REF], thus questioning the underlying assumption of sparsity often used to estimate brain network connectivity. On the other hand, for the problem of finding differences between networks in two populations, sparsity may be a valid assumption. It is well known that anatomical brain connections tend to closely follow functional ones [START_REF] Honey | Predicting human resting-state functional connectivity from structural connectivity[END_REF]. Since anatomical networks do not differ drastically we can surmise that two brain networks should not differ much even in the presence of pathologies. The statistical method we present here leverages sparsity in the difference of two networks, to yield well-behaved estimation and hypothesis testing in the low-sample regime. Most closely related to our work, [35,[START_REF] Fazayeli | Generalized direct change estimation in ising model structure[END_REF] recently consider a different approach to estimating difference networks, but does not consider assigning significance to the detection of edges.

Debiased Difference Estimation

In many applications one may be interested in learning multiple linear models from data that share many parameters. Situations such as this arise often in neuroimaging and bioinformatics applications. We can often improve the learning procedure of such models by incorporating fused penalties that penalize the • 1 norm of the parameter differences or • 1,2 which encourages groups of parameters to shrink together. These methods have been shown to substantially improve the learning of the joint models. However, the differences between model parameters, which can have a high sample complexity when there are few of them, are often pointed out only in passing [START_REF] Chen | Smoothing proximal gradient method for general structured sparse learning[END_REF][START_REF] Danaher | The joint graphical lasso for inverse covariance estimation across multiple classes[END_REF][START_REF] Honorio | Multi-task learning of Gaussian graphical models[END_REF]. On the other hand, in many situations we might be interested in actually understanding and identifying the differences between elements of the support. For example when considering brain networks of patients suffering from a pathology and healthy control subjects, the difference in brain connectivity may be of great interest. Here we focus specifically on accurately identifying differences with significance.

We consider the case of two tasks (e.g. two groups of subjects), but the analysis can be easily extended to general multi-task settings. Consider the problem setting of data matrices X 1 and X 2 , which are n 1 × p and n 2 × p, respectively. We model them as producing outputs Y 1 and Y 2 , corrupted by diagonal gaussian noise 1 and 2 as follows

Y 1 = X 1 β 1 + 1 , Y 2 = X 2 β 2 + 2 (3)
Let S 1 and S 2 index the elements of the support of β 1 and β 2 , respectively. Furthermore the support of β 1 -β 2 is indexed by S d and finally the union of S 1 and S 2 is denoted S a . Using a squared loss estimator producing independent estimates β1 , β2 we can obtain a difference estimate βd = β1 -β2 . In general if S d is very small relative to S a then we will have a difficult time to identify the support S d . This can be seen if we consider each of the individual components of the prediction errors. The larger the true support S a the more it will drown out the subset which corresponds to the difference support. This can be true even if one uses 1 regularizers over the parameter vectors. Consequently, one cannot rely on the straightforward strategy of learning two independent estimates and taking their difference. The problem is particularly pronounced in the common setting where one group has fewer samples than the other. Thus here we consider the setting where n 1 > n 2 and possibly n 1 n 2 .

Let β1 and β2 be regularized least squares estimates. In our problem setting we wish to obtain confidence intervals on debiased versions of the difference βd = β1 -β2 in a high-dimensional setting (in the sense that n 2 < p), we aim to leverage assumptions about the form of the true β d , primarily that it is sparse, while the independent β1 and β2 are weakly sparse or not sparse. We consider a general case of a joint regularized least squares estimation of β1 and β2 min β1,β2

1 n 1 Y 1 -X 1 β 1 2 + 1 n 2 Y 2 -X 2 β 2 2 + R(β 1 , β 2 ) (4)
We note that the differentiating and using the KKT conditions gives

kλ = k1 k2 = 1 n1 X T 1 (Y -X 1 β 1 ) 1 n2 X T 2 (Y -X 2 β 2 ) (5)
where kλ is the (sub)gradient of R(β 1 , β 2 ). Substituting Equation ( 3) we can now write

Σ1 ( β1 -β 1 ) + k1 = 1 n 1 X T 1 1 and Σ2 ( β2 -β 2 ) + k2 = 1 n 2 X T 2 2 (6)
We would like to solve for the difference β1 -β2 but the covariance matrices may not be invertible.

We introduce matrices M 1 and M 2 , which will allow us to isolate the relevant term. We will see that in addition these matrices will allow us to decouple the bias and variance of the estimators.

M 1 Σ1 ( β1 -β 1 ) + M 1 k1 = 1 n 1 M 1 X T 1 1 and M 2 Σ2 ( β2 -β 2 ) + M 2 k2 = 1 n 2 M 2 X T 2 2 (7)
subtracting these and rearranging we can now isolate the difference estimator plus a term we add back controlled by M 1 and M 2

( β1 -β2 ) -(β 1 -β 2 ) + M 1 k1 -M 2 k2 = 1 n 1 M 1 X T 1 1 - 1 n 2 M 2 X T 2 2 -∆ (8)
∆ = (M 1 Σ1 -I)( β1 -β 1 ) -(M 2 Σ2 -I)( β2 -β 2 ) (9)
Denoting

β d := β 1 -β 2 and β a := β 1 + β 2 , we can reformulate ∆: ∆ = (M 1 Σ1 -I + M 2 Σ2 -I) 2 ( βd -β d ) + (M 1 Σ1 -M 2 Σ2 ) 2 ( βa -β a ) (10)
Here, ∆ will control the bias of our estimator. Additionally, we want to minimize its variance,

1 n 1 M 1 Σ1 M 1 σ2 1 + 1 n 2 M 2 Σ2 M 2 σ2 2 . ( 11
)
We can now overcome the limitations of simple bias variance trade-off by using an appropriate regularizer coupled with an assumption on the underlying signal β 1 and β 2 . This will in turn make ∆ asymptotically vanish while maximizing the variance.

Since we are interested in pointwise estimates, we can focus on bounding the infinity norm of ∆.

∆ ∞ ≤ 1 2 M 1 Σ1 + M 2 Σ2 -2I ∞ µ1 βd -β d 1 l d + 1 2 M 1 Σ1 -M 2 Σ2 ∞ µ2 βa -β a 1 la (12)
We can control the maximum bias by selecting M 1 and M 2 appropriately. If we use an appropriate regularizer coupled with sparsity assumptions we can bound the terms l a and l d and use this knowledge to appropriately select M 1 and M 2 such that the bias becomes neglibile. If we had only the independent parameter sparsity assumption we can apply the results of the debiased lasso and estimate M 1 and M 2 independently as in [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF]. In the case of interest where β 1 and β 2 share many weights we can do better by taking this as an assumption and applying a sparsity regularization on the difference by adding the term λ 2 β 1 -β 2 1 . Comparing the decoupled penalty to the fused penalty proposed we see that l d would decrease at a given sample size. We now show how to jointly estimate M 1 and M 2 so that ∆ ∞ becomes negligible for a given n, p and sparsity assumption.

Debiasing the Multi-Task Fused Lasso

Motivated by the inductive hypothesis from neuroscience described above we introduce a consistent low-variance estimator, the debiased multi-task fused lasso. We propose to use the following regularizer R(

β 1 , β 2 ) = λ 1 β 1 1 + λ 1 β 2 1 + λ 2 β 1 -β 2 1
. This penalty has been referred to in some literature as the multi-task fused lasso [START_REF] Chen | Smoothing proximal gradient method for general structured sparse learning[END_REF]. We propose to then debias this estimate as shown in [START_REF] Martino | The autism brain imaging data exchange: Towards a large-scale evaluation of the intrinsic brain architecture in autism[END_REF]. We estimate the M 1 and M 2 matrices by solving the following QP for each row m 1 and m 2 of the matrices M 1 and M 2 .

min m1,m2 1 n 1 m T 1 Σ1 m 1 + 1 n 2 m T 2 Σ2 m 2 (13)
s.t. M 1 Σ1 + M 2 Σ2 -2I ∞ ≤ µ 1 , M 1 Σ1 -M 2 Σ2 ∞ ≤ µ 2
This directly minimizes the variance, while bounding the bias in the constraint. We now show how to set the bounds:

Proposition 1. Take λ 1 > 2 log p n2 and λ 2 = O(λ 1 ). Denote s d the difference sparsity, s 1,2 the parameter sparsity |S 1 | + |S 2 |, c > 1,a > 1, and 0 < m 1.
When the compatibility condition [START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF][START_REF] Ganguly | Local neighborhood fusion in locally constant Gaussian graphical models[END_REF] holds the following bounds gives l a u 2 = o(1) and l d u 1 = o(1) and thus ∆ ∞ = o(1) with high probability.

µ 1 ≤ 1 cλ 2 s d n m 2 and µ 2 ≤ 1 a(λ 1 s 1,2 + λ 2 s d )n m 2 ( 14
)
The proof is given in the supplementary material. Using the prescribed M s obtained with (13) and 14 we obtain an unbiased estimator given by [START_REF] Martino | The autism brain imaging data exchange: Towards a large-scale evaluation of the intrinsic brain architecture in autism[END_REF] with variance (11)

GGM Difference Structure Discovery with Significance

The debiased lasso and the debiased multi-task fused lasso, proposed in the previous section, can be used to learn the structure of a difference of Gaussian graphical models and to provide significance results on the presence of edges within the difference graph. We refer to these two procedures as Difference of Neighborhoods Debiased Lasso Selection and Difference of Neighborhoods Debiased Fused Lasso Selection.

We recall that the conditional independence properties of a GGM are given by the zeros of the precision matrix and these zeros correspond to the zeros of regression parameters when regressing one variable on all the other. By obtaining a debiased lasso estimate for each node in the graph [START_REF] Waldorp | Testing for graph differences using the desparsified lasso in high-dimensional data[END_REF] notes this leads to a sparse unbiased precision matrix estimate with a known asymptotic distribution. Subtracting these estimates for two different datasets gives us a difference estimate whose zeros correspond to no difference of graph edges in two GGMs. We can similarly use the debiased multitask fused lasso described above and the joint debiasing procedure to obtain a test statistic for the difference of networks. We now formalize this procedure.

Notation Given GGMs j = 1, 2. Let X j denote the random variable in R p associated with GGM j. We denote X j,v the random variable associated with a node, v of the GGM and X j,v c all other nodes in the graph. We denote βj,v the lasso or multi-task fused lasso estimate of X j,v c onto X j,v , then βj,dL,v is the debiased version of βj,v . Finally let β j,v denote the unknown regression, X j,v = X j,v c β j,v + j where j ∼ N(0, σ j I). Define β i D,v = βi 1,dL,v -βi 2,dL,v the test statistic associated with the edge v, i in the difference of GGMs j = 1, 2.

Algorithm 1 Difference Network Selection with

Neighborhood Debiased Lasso V = {1, ..., P } NxP Data Matrices, X1 and X2 Px(P-1) Output Matrix B of test statistics for v ∈ V do Estimate unbiased σ1, σ2 from X1,v, X2,v for j ∈ {1, 2} do βj ← SolveLasso(Xj,vc , Xj,v)

Mj ← M Estimator(Xj,vc ) βj,U ← βj+MjX T j,v c (Xj,v -Xj,vc βj) end for σ 2 d ← diag( σ2 1 n 1 M T 1 Σ1M1 + σ2 2 n 2 M T 2 Σ2M2) for j ∈ v c do Bv,j = (β1,U,j -β2,U,j)/ σ
for j ∈ {1, 2} do βj,U ← βj+MjX T j,v c (Xj,v -Xj,vc βj) end for σ 2 d ← diag( σ2 1 n 1 M T 1 Σ1M1 + σ2 2 n 2 M T 2 Σ2M2) for j ∈ v c do Bv,j = (β1,U,j -β2,U,j)/ σ 2 d,j
end for end for Proposition 2. Given the βi D,v , M 1 and M 2 computed as in [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF] for the debiased lasso or as in Section 3.1 for the debiased multi-task fused lasso. When the respective assumptions of these estimators are satisfied the following holds w.h.p.

βi D,v -β i D,v = W + o(1) where W ∼ N(0, [σ 2 1 M 1 Σ1 M T 1 + σ 2 2 M 2 Σ2 M T 2 ] i,i ) (15)
This follows directly from the asymptotic consistency of each individual βi j,dL,v for the debiased lasso and multi-task fused lasso.

We can now define the the null hypothesis of interest as H 0 : Θ 1,(i,j) = Θ 2,(i,j) . Obtaining a test statistic for each element β i D,v allows us to perform hypothesis testing on individual edges, all the edges, or groups of edges (controlling for the FWER). We summarize the Neighbourhood Debiased Lasso Selection process in Algorithm 1 and the Neighbourhood Debiased Multi-Task Fused Lasso Selection in Algorithm 2 which can be used to obtain a matrix of all the relevant test statistics.

Experiments

Simulations

We generate synthetic data based on two Gaussian graphical models with 75 vertices. Each of the individual graphs have a sparsity of 19% and their difference sparsity is 3%. We construct the models by taking two identical precision matrices and randomly removing some edges from both. We generate synthetic data using both precision matrices. We use n 1 = 800 samples for the first dataset and vary the second dataset n 2 = 20, 30, ...150.

We perform a regression using the debiased lasso and the debiased multi-task fused lasso on each node of the graphs. As an extra baseline we consider the projected ridge method from the R package "hdi" [START_REF] Dezeure | High-dimensional inference: Confidence intervals, p-values and R-software hdi[END_REF]. We use the debiased lasso of [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF] where we set λ = kσ log p/n. We select c by 3fold cross validation k = {0.1, ..100} and M as prescribed in [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF] which we obtain by solving a quadratic program. σ is an unbiased estimator of the noise variance. For the debiased lasso we let both λ 1 = k 1 σ2 log p/n 2 and λ 2 = k 2 σ2 log p/n 2 , and select based on 3-fold cross-validation from the same range as k. M 1 and M 2 are obtained as in Equation ( 13) with the bounds (14) being set with c = a = 2, s d = 2, s 1,2 = 15, m = 0.01, and the cross validated λ 1 and λ 2 . In both debiased lasso and fused multi-task lasso cases we utilize the Mosek QP solver package to obtain M . For the projected ridge method we use the hdi package to obtain two estimates of β 1 and β 2 along with their upper bounded biases which are then used to obtain p-values for the difference.

We report the false positive rate, the power, the coverage and interval length as per [START_REF] Van De Geer | On asymptotically optimal confidence regions and tests for high-dimensional models[END_REF] for the difference of graphs. In these experiments we aggregate statistics to demonstrate power of the test statistic, as such we consider each edge as a separate test and do not perform corrections. Table 1 gives the numerical results for n 2 = 60: the power and coverage is substantially better for the debiased fused multi-task lasso, while at the same time the confidence interval smaller. Fused Lasso, and Projected Ridge Regression for edge selection in difference of GGM. The significance level is 5%, n1 = 800 and n2 = 60. All methods have false positive below the significance level and the debiased fused lasso dominates in terms of power. The coverage of the difference support and non-difference support is also best for the debiased fused lasso, which simultaneously has smaller confidence intervals on average.

Figure 1 shows the power of the test for different values of n 2 . The fusedlasso outperforms the other methods substantially. Projected ridge regression is particularly weak, in this scenario, as it uses a worst case p-value obtained using an estimate of an upper bound on the bias [START_REF] Dezeure | High-dimensional inference: Confidence intervals, p-values and R-software hdi[END_REF].

Autism Dataset

Correlations in brain activity measured via fMRI reveal functional interactions between remote brain regions [START_REF] Lindquist | The statistical analysis of fMRI data[END_REF]. In population analysis, they are used to measure how connectivity varies between different groups. Such analysis of brain function is particularly important in psychiatric diseases, that have no known anatomical support: the brain functions in a pathological aspect, but nothing abnormal is clearly visible in the brain tissues. Autism spectrum disorder is a typical example of such ill-understood psychiatric disease. Resting-state fMRI is accumulated in an effort to shed light on this diseases mechanisms: comparing the connectivity of autism patients versus control subjects. The ABIDE (Autism Brain Imaging Data Exchange) dataset [START_REF] Martino | The autism brain imaging data exchange: Towards a large-scale evaluation of the intrinsic brain architecture in autism[END_REF] gathers rest-fMRI from 1 112 subjects across, with 539 individuals suffering from autism spectrum disorder and 573 typical controls. We use the preprocessed and curated data 1 .

In a connectome analysis [START_REF] Varoquaux | Learning and comparing functional connectomes across subjects[END_REF][START_REF] Richiardi | Decoding brain states from fMRI connectivity graphs[END_REF], each subject is described by a GGM measuring functional connectivity between a set of regions. We build a connectome from brain regions of interest based on a multi-subject atlas2 of 39 functional regions derived from resting-state fMRI [START_REF] Varoquaux | Multi-subject dictionary learning to segment an atlas of brain spontaneous activity[END_REF] (see. Fig. 4).

We are interested in determining edge differences between the autism group and the control group. We use this data to show how our parametric hypothesis test can be used to determine differences in brain networks. Since no ground truth exists for this problem, we use permutation testing to evaluate the statistical procedures [START_REF] Nichols | Nonparametric permutation tests for functional neuroimaging: A primer with examples[END_REF][START_REF] Mota | Randomized parcellation based inference[END_REF]. Here we permute the two conditions (e.g. autism and control group) to compute a p-value and compare it to our test statistics. This provides us with a finite sample strict control on the error rate: a non-parametric validation of our parametric test.

For our experiments we take 2000 randomly chosen volumes from the control group subjects and 100 volumes from the autism group subjects. We perform permutation testing using the de-biased lasso, de-biased multi-task fused lasso, and projected ridge regression. Parameters for the de-biased fused lasso are chosen as in the previous section. For the de-biased lasso we use the exact settings for λ and constraints on M provided in the experimental section of [START_REF] Javanmard | Confidence intervals and hypothesis testing for high-dimensional regression[END_REF]. Projected ridge regression is evaluated as in the previous section.

Figure 2 shows a comparison of three parametric approaches versus their analogue obtained with a permutation test. The chart plots the permutation p-values of each entry in the 38 × 39 B matrix against the expected parametric p-value. For all the methods the points are above the line indicating the tests are not breaching the expected false positive rates. However the de-biased lasso and ridge projecting are very conservative and lead to few detections. The de-biased multi-task fused lasso yields far more detections on the same dataset, within the expected false positive rate or near it.

We now analyse the reproducibility of the results by repeatedly sampling 100 subsets of the data (with the same proportions n 1 = 2000 and n 2 = 100), obtaining the matrix of test statistics, selecting edges that fall below the 5% significance level. Figure 3 shows how often edges are selected multiple times across subsamples. We report results with a threshold on uncorrected p-values as the lasso procedure selects no edges with multiple comparison correction (supplementary materials give FDR-corrected results for the de-biased fused multi-task lasso selection). Figure 5 shows a connectome of the edges frequently selected by the de-biased fused multi-task lasso (with FDR correction).

Conclusions

We have shown how to characterize the distribution of differences of sparse estimators and how to use this distribution for confidence intervals and p-values on GGM network differences. For this purpose, we have introduced the de-biased multi-task fused lasso. We have demonstrated on synthetic and real data that this approach can provide accurate p-values and a sizable increase of statistical power compared to standard procedures. The settings match those of population analysis for functional brain connectivity, and the gain in statistical power is direly needed to tackle the low sample sizes [START_REF] Button | Power failure: Why small sample size undermines the reliability of neuroscience[END_REF].

Future work calls for expanding the analysis to cases with more than two groups as well as considering a 1,2 penalty sometimes used at the group level [START_REF] Varoquaux | Brain covariance selection: Better individual functional connectivity models using population prior[END_REF]. Additionally the squared loss objective optimizes excessively the prediction and could be modified to lower further the sample complexity in terms of parameter estimation.

Proof. We first consider the bound associated with la λ1 βa -βa 1 ≤ λ1 βS

In the domain of interest n1 n2 if we select λ2 = O(λ1) we can see the relevant terms related to the parameter support become small with respect to terms with S1,2. Thus the error on the difference should dominate. In this region we can have 3λ1 βS 1,2 -βS 1,2 1/2 -λ1 βS c 1,2 1/2 ≤ cλ2 βd,S -β d,S 1 where c > 0.

B Additional Experimental Details

We show the corrected reproducibility results in Figure 6. For multiple testing correction in our experiments We use the Benjamin-Hochberg FDR procedure. Reproducibility across subsamples lasso fusedlasso Figure 6: Reproducibility of results from sub-sampling using FDR of 5% Reproducibility of results from subsampling, debiased lasso does not produce any significant edge differences that correspond to a 5% error rate

Figure 1 :

 1 Figure 1: Power of the test for different number of samples in the second simulation, with n1 = 800. The debiased fused lasso has highest statistical power.

Figure 2 :

 2 Figure 2: Permutation testing comparing debiased fused lasso, debiased lasso, and projected ridge regression on the ABIDE dataset. The chart plots the permutation p-values of each method on each possible edge against the expected parametric p-value. The debiased lasso and ridge projection are very conservative and lead to few detections. The fused lasso yields far more detections on the same dataset, almost all within the expected false positive rate.

  Fraction of connections occuring at least t timesReproducibility across subsamples lasso fusedlasso

Figure 3 :

 3 Figure 3: Reproducibility of results from sub-sampling using uncorrected error rate. The fused lasso is much more likely to detect edges and produce stable results. Using corrected p-values no detections are made by lasso (figure in supplementary material).

Figure 4 :Figure 5 :

 45 Figure 4: Outlines of the regions of the MSDL atlas.L R L R fusedlasso_subsamples_100

λ2 βd -β d 1 ≤ 2 ( 33 )≤ c 2 λ 2 2 s d φ 2 min+Thus βd -β d 1 ≤ c 2 λ 2 s d φ 2 min

 1233212 2λ2 βd,S -β d,S 1 -XN ( β -β) βd,S -β d,S 2 -XN ( β -β) 2Invoking the compatibility assumption[START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF] ≤ 2cλ2 √ s d XN ( β -β) 2/φmin -XN ( β -β) XN ( β -β) 2 2 -XN ( β -β)and use of the bound prescribed gives l d u1 = o(1).

  Fraction of connections occuring at least t times

Table 1 :

 1 Comparison of Debiased Lasso, Debiased

	Method	FP	TP(Power) Cov S Cov S c d	len S len S c d
	Deb. Lasso	3.7%	80.6%	96.2%	92%	2.199	2.195
	Deb. Fused Lasso 0.0%	93.3%	100%	98.6% 2.191	2.041
	Ridge Projection	0.0%	18.6%	100%	100%	5.544	5.544

  1,2 -βS 1,2 1 + λ1 βS c 1,2 1 ≤ 4λ1 βS 1,2 -βS 1,2 1 + 2λ2 βS d -βS d 1 -2 XN ( β -β) 2 βS 1,2 -βS 1,2 2 + 2λ2 √ s d βS d -βS d 2The bound u2 now follows by inverting the expression shown and adding a factor of n m 2 where m 1.Now we consider the bound for l d .λ2 βd -β d 1 = λ2 βd,S -β d,S 1 + λ2 βd,S c 1(29)≤ 2λ2 βd,S -β d,S 1 + 3λ1 βS 1,2 -βS 1,2 1/2

						2	(25)
	≤ 4λ1	√ s1,2 -2 XN ( β -β) 2 2	(26)
	Invoking the compatibility assumption [1, 16, 11] with compatibility constant φmin
	≤	4λ1 φmin √ s1,2	XN ( β -β) 2 +	2λ2 φmin √ s d	XN ( β -β) 2
				-2 XN ( β -β) 2 2	(27)
				≤	4λ 2 1 s1,2 φ 2 min	+	2λ 2 2 s d min φ 2	(28)
						(30)
			-XN ( β -β) 2 2 -λ1 βS c 1,2 1/2

http://preprocessed-connectomes-project.github.io/abide/

https://team.inria.fr/parietal/research/spatial_patterns/ spatial-patterns-in-resting-state/

Acknowledgements

This work is partially funded by Internal Funds KU Leuven, ERC Grant 259112, FP7-MC-CIG 334380, and DIGITEO 2013-0788D -SOPRANO, and ANR-11-BINF-0004 NiConnect. We thank Wacha Bounliphone and Puneet Dokania for helpful discussions

A Analysis of Debiased Multi-Task Fused Lasso

The following analysis is used to show the conditions under which the debias multi-task fused lasso achieves a negligible bias.

This follows from the fact that β is the minimizer of the fused lasso objective.

The term, T XN ( β -β), commonly known as the empirical process term [START_REF] Bühlmann | Statistics for High-Dimensional Data[END_REF] can be bound as follows:

Where we utilize holder's inequality in the last line. We define the random event F for which the following holds:

Lemma 2. Suppose Σj,j = 1 for both X1 and X2 then we have for all t > 0 and n1 > n2

Proof. This follows directly from the [1, Lemma 6.2] and taking n1 > n2.

This allows us to get rid of the empirical process term on F, with an appropriate choice of λ1.

Given a set, S, denote βS the vector of equal size to β but all elements not in S set to zero. We can now show the following Lemma 3. We have on F with λ1 ≥ 2λ0

Proof. Following [1, Lemma 6.3] we start with the basic inequality on F. Which gives

Since we assume the truth is in fact sparse,

Furthermore,

Substituting ( 22), [START_REF] Mohan | Structured learning of Gaussian graphical models[END_REF], and ( 21) into [START_REF] Lockhart | A significance test for the lasso[END_REF] and rearranging completes the proof.

From the lemma above we can now justify the bounds in ( 14