Multilevel Particle Filters for Partially Observed McKean-Vlasov Stochastic Differential Equations
Elsiddig Awadelkarim Ajay Jasra






verbose,tmargin=3cm,bmargin=3cm,lmargin=2.75cm,rmargin=2.75cm





































































ifundefineddefinecolorifundefineddefinecolor
 ifundefineddefinecolor
 









assAssumption[section]










theoremTheorem[section]
lemLemma[section]
remRemark[section]
propProposition[section]
corCorollary[section]
hypA
hypAhypA
  
  (AhypA)
hypB
hypBhypB
 
  (HhypB)
hypD
hypDhypD
 
  (DhypD)
hypW
hypWhypW
  
  (WhypW)




















































































1016

 








Multilevel Particle Filters for Partially Observed McKean-Vlasov Stochastic Differential Equations




BY  ELSIDDIG AWADELKARIM ^1 & AJAY JASRA^2

^1Applied Mathematics and Computational Science Program,  Computer, Electrical and Mathematical Sciences and Engineering Division, King Abdullah University of Science and Technology, Thuwal, 23955-6900, KSA.

^2School of Data Science, The Chinese University of Hong Kong, Shenzhen, CN.
E-Mail:  





In this paper we consider the filtering problem associated to partially observed McKean-Vlasov stochastic differential equations (SDEs). The model consists of data that are observed at regular and discrete times and the objective is to compute the conditional expectation of (functionals) of the solutions of the SDE at the current time.
This problem, even the ordinary SDE case is challenging and requires numerical approximations. Based upon the ideas in <cit.> we develop a new particle filter (PF) and multilevel particle filter (MLPF) to approximate the afore-mentioned expectations. We prove under assumptions that, for ϵ>0, to obtain a mean square error of 𝒪(ϵ^2) the PF has a cost per-observation time of 𝒪(ϵ^-5) and the MLPF costs 𝒪(ϵ^-4) (best case) or 𝒪(ϵ^-4log(ϵ)^2) (worst case). Our theoretical results are supported by numerical experiments.


Key words: McKean-Vlasov SDEs,  Importance Sampling,  Multilevel Particle Filters.





§ INTRODUCTION


We consider the stochastic differential equation (SDE) with X_0=x_0∈ℝ^d fixed:

    dX_t = a(X_t,ξ_1(X_t,μ_t))dt + b(X_t,ξ_2(X_t,μ_t))dW_t

where for j∈{1,2}

    ξ_j(X_t,μ_t)  =  ∫_ℝ^dξ_j(X_t,x)μ_t(dx)

{W_t}_t≥ 0 is a standard d-dimensional Brownian motion,  ξ_j:ℝ^2d→ℝ, a:ℝ^d×ℝ→ℝ^d, b:ℝ^d×ℝ→ℝ^d× d and μ_t is the law of X_t.  This is the collection McKean-Vlasov SDEs
and finds a wide variety of real applications such as biology and finance <cit.>, in addition they have been shown to be equivalent representation of solutions of classes of nonlinear partial differential equations <cit.>,  and one can expect even more applications in the future. 
The exact problem that will be addressed in this article is where the process defined in (<ref>) is only observed with noise at regular and discrete times. The objective is then to compute the conditional expectation of X_t (or functionals thereof) given all the observations recursively at every observation time, which is a filtering problem; see e.g. <cit.> for an introduction.

The filtering of continuous-time processes is very challenging. First, even if the transition densities of the process exist and are analytically tractable,  filtering must normally be perfomed by using numerical methods of which we focus upon particle filters (PF); see e.g. <cit.>. Second,  the afore-mentioned transition densities are seldom available and hence one normally has to resort to a time discretization of the process; see for instance <cit.>. In the case of regular diffusions (i.e. where the coefficients do no depend on μ_t), there is already a substantial literature on using PFs for solving the filtering problem; see e.g. <cit.>.  In addition to these methods PFs have been extended to leverage on the popular multilevel Monte Carlo (MLMC) method <cit.>.
These multilevel PFs (MLPFs) <cit.> can reduce the cost to compute the filtering expectation for a given mean square error (MSE).

To the best of our knowledge,  there is little work in the context of filtering McKean-Valasov SDEs. However, when
trying to compute expectations with respect to (w.r.t.) the laws of such processes, several methods have been proposed. This includes the initial work of <cit.> that relies on time discretization and simulation coupled with solving backward Kolmogorov equations; this can be prohibitive if d is even moderate.  Several recent extensions based upon using the afore-mentioned simulation and designing change of measures for importance sampling have been considered in <cit.> and then combining with MLMC <cit.> or multi-index MC (see <cit.>) <cit.>.  An iterative method, combined with MLMC was also proposed in <cit.>.  Despite these high quality and interesting contributions, none of these methods are designed for the filtering problem and would likely have several deficiencies.  For instance,  if the change of measure methods would be adopted it is likely that they would suffer from the well-known weight-degeneracy problem (e.g. <cit.>). That is,  the estimates of the filter would suffer an exponential increase of variance w.r.t. to observation times. The approach of <cit.> is iterative and thus is not appropriate for the filtering problem as one seeks fast (non-iterative) methods to assimilate the data.

In this paper we develop new PF and MLPF methodology to estimate the filter, associated to partially observed Mckean-Vlasov SDEs, recursively in time. This requires some care as, unlike the case of filtering ordinary diffusions, even under a time-discretization, one cannot simulate the SDE. This is because the law of the discretized process is not known and must be estimated. We show how the idea adopted in <cit.> can be leveraged. That is, to use the particle system of <cit.> to be fed into a `classical' PF and this resulting algorithm can be used to recursively estimate the filter. We then show how this idea can be extended to the multilevel context by combining with ideas from <cit.>; see <cit.> for a review and <cit.> for related methodology.
We prove, under assumptions, that for, ϵ>0 to achieve a MSE, associated to estimating the filter, of 𝒪(ϵ^2) per-observation time:


  * the cost for the PF is 𝒪(ϵ^-5)

  * the cost for the MLPF is 𝒪(ϵ^-4) (best case) or 𝒪(ϵ^-4log(ϵ)^2) (worst case).

The best case corresponds to the case where the diffusion coefficient b in (<ref>) is constant and the other case is for the non-constant case - of course only when the SDE follows our assumptions.
The costs above are higher than one might expect for ordinary Monte Carlo based methods (of which PFs and MLPFs are). These increased costs are primarily associated to having to approximate the laws μ_t (and their time-discretized versions) in the evolution of the SDE.  There are possible ways to deal with this, but it is not investigated in this paper. The MSE to cost results are verified in several numerical examples.
























This article is structured as follows. In Section <ref> we give our approach to filtering partially observed McKean-Vlasov SDEs.  In Section <ref> our theoretical results are presented.  Section <ref> details our numerical results associated to our methodology.  Our proofs can be found in the Appendix.



§ APPROACH




 §.§ State-Space Model


We denote by P_μ_t-1,t(x_t-1,dx_t) the conditional law of X_t (as given in (<ref>)) given ℱ_t-1 (the natual filtration of the process), for t≥ 1; that is, the transition kernel over unit time. We consider a discrete time observation
process Y_1,Y_2,…,  Y_t∈𝖸, that are assumed, for notational convenience, to be observed at unit times. Conditional on the position X_t, t∈ℕ of (<ref>), the random variable Y_t is assumed to be independent of all other random variables, with a bounded and positive probability density G(x_t,y_t).

Let φ:ℝ^d→ℝ be a bounded and measurable function (write such a collection of functions ℬ_b(ℝ^d)), then we define the filtering expectation for t∈ℕ as:

    π_t(φ) := ∫_ℝ^dtφ(x_t){∏_p=1^t G(x_p,y_p)}∏_p=1^tP_μ_p-1,p(x_p-1,dx_p)/∫_ℝ^dt{∏_p=1^t G(x_p,y_p)}∏_p=1^tP_μ_p-1,p(x_p-1,dx_p)

where μ_0=δ_{x_0} (Dirac measure on the set {x_0}). We remark that one does not need φ
and each of the G(·,y_t) to be bounded, but it will simplify the resulting exposition to do so.

In most cases of practical interest, μ_t and the dynamics P_μ_t-1,t(x_t-1,dx_t) are difficult to work with.  For instance the transition kernel cannot be simulated in many problems. As a first step to approximate the filter π_t, we introduce a time-discretization over a regular grid of spacing Δ_l=2^-l, l∈ℕ_0. The significance of subscript l is that it will denote the level of discretization - as l grows so will the accuracy of the time discretization. We will use the Euler-Maruyama method associated to (<ref>) and denote the law at any time t∈{0,Δ_l,2Δ_l,…} as μ_t^l. That is, we now consider the approximation for k∈ℕ_0:

    X_kΔ_l    =     X_(k-1)Δ_l + a(X_(k-1)Δ_l,ξ_1(X_(k-1)Δ_l,μ_(k-1)Δ_l^l)) +
          b(X_(k-1)Δ_l,ξ_2(X_(k-1)Δ_l,μ_(k-1)Δ_l^l))[W_kΔ_l - W_(k-1)Δ_l]

where X_0=x_0 and μ_0^l=δ_{x_0}.  Associated to (<ref>), we denote by P_μ_t-1,t^l(x_t-1,dx_t) the conditional law of X_t, t∈ℕ,  given ℱ_t-1 for t≥ 1; that is, the transition kernel over unit time induced by (<ref>).  It should be remarked that in many cases, 
(<ref>) cannot be simuated exactly as the expectations associated to μ_(k-1)Δ_l^l cannot be computed even if one knows μ_(k-1)Δ_l^l, which is again unlikely.

Our objective now is to compute an approximation of the time discretized filtering expectation for (t,l,φ)∈ℕ×ℕ_0×ℬ_b(ℝ^d):

    π_t^l(φ) := ∫_ℝ^dtφ(x_t){∏_p=1^t G(x_p,y_p)}∏_p=1^tP_μ_p-1^l,p^l(x_p-1,dx_p)/∫_ℝ^dt{∏_p=1^t G(x_p,y_p)}∏_p=1^tP_μ_p-1^l,p^l(x_p-1,dx_p).





 §.§ Approximating the Law


We now consider a method that will be used to provide a Monte Carlo based approximation of the law
μ_t^l. The approach we present is a simple discretized method in Algorithm <ref> from <cit.>.   In Algorithm <ref>, the notation 𝒩_d(κ,Σ) denotes the d-dimensional Gaussian distribution with mean κ and covariance matrix Σ. I_d is the d× d
identity matrix and ind∼ denotes independently distributed as.

Algorithm <ref> can be used to approximate expectations w.r.t. μ_t^l
and indeed on the grid in-between time t-1 and t.   Algorithm <ref> is given in the form that we need it later on in the article.





 §.§ Particle Filter


An approach for the recursive approximation of π_t^l(φ) is
 presented in Algorithm <ref> and is run to a terminal time T - this need not be the case.
At this stage, several remarks are of interest:


  * The cost of Algorithm <ref> is 𝒪(Δ_l^-1M(M+N)) per observation time (or time step).

  * In step 2. it is critical that this process be independent of the particles that are generated in the particle filter.  The empirical laws produced here are generated independently of anything that is done later on.
If one used the empirical laws that were produced by the samples from the particle filter, then the approximation would be biased even if N,M→∞ and l→∞. This is because these particles used in the particle filter have been designed to approximate expectations w.r.t. π_t and not the law of the McKean-Vlasov SDE itself. Thus, here the motivation of a driving approximation of the SDE marginal law, as was the case in <cit.>,  is one of necessity. 

  * The variance w.r.t. the time parameter should be controlled and can be uniform in time.
This may not be the case for the methods in <cit.>.

  * Improved discretization methods can be used, when available.

  * We will prove, later on, that as N,M→∞ we have that π_t^l,N,M(φ)→_ℙπ_t^l(φ), where →_ℙ denotes convergence in probability and φ is in a particular class of real-valued functions to be stated.








 §.§ Multilevel Particle Filter


We now show how the PF can be extended into the MLMC framework,  along the lines of the work in <cit.>.
The purpose of this is to reduce the computing cost in terms of an MSE, which is to be defined in Section <ref>.  The basic idea is to approximate the identity for (L,t,φ)∈ℕ^2×ℬ_b(ℝ^d):

    π_t^L(φ) = π_t^0(φ) + ∑_l=1^L {π_t^l(φ)-π_t^l-1(φ)}.

The terms on the right hand side (R.H.S.) of the above equation are then approximated as follows. The term
π_t^0(φ) can be dealt with using the PF of Algorithm <ref>. Then independently of using the PF and for each summand, we run what is called a coupled PF (CPF). This is a type of PF which allows one to approximate 
π_t^l(φ)-π_t^l-1(φ) by producing a coupling of the sampling step (Algorithm <ref> step 3.)
and a coupling of the resampling step (Algorithm <ref> step 4.).  
In Algorithm <ref> we also use Algorithm <ref> which needs to be coupled as well.
This coupling and also for the sampling step of a PF is achieved using the well-known synchronous coupling for diffusions (and the associated Euler approximation) and is given in our coupling of Algorithm <ref>, which can be found in Algorithm <ref>.
The coupling of the resampling step relies on simulating a maximal coupling and is given in
Algorithm <ref>.  As noted in <cit.>, there is nothing optimal about using this approach in a CPF,  but it is rather popular in the literature and is of linear complexity in terms of the cardinality of the state-space.

Given all the ideas already stated,  we are now in a position to give the CPF algorithm which is Algorithm <ref>.  The cost of the CPF is 𝒪(Δ_l^-1M(M+N)) per observation time.
We can also describe the MLPF method that we will use for the approximation of π_t^L(φ).
We assume that we have L∈ℕ given and for each level l, the number of samples N_l used at each level of the application of Algorithm <ref>, with the number of samples used for the call of 
Algorithm <ref> within Algorithm <ref> equal to M_l also; we shall describe how to choose these numbers later on.  Then the procedure would be as follows:


  * Run Algorithm <ref> with N_0 samples and M_0 samples for each call of Algorithm <ref>.

  * Independently for each l∈{1,…,L}, run Algorithm <ref> with N_l samples and M_l samples for each call of Algorithm <ref>.

Then our MLPF estimator of π_t^L(φ) is for (L,t,φ)∈ℕ^2×ℬ_b(ℝ^d):

    π_t^L(φ) = π_t^0,N_0,M_0(φ) + ∑_l=1^L {π_t^l,N_l,M_l(φ) - π_t^l-1,N_l,M_l(φ)}

where π_t^0,N_0,M_0(φ) has been computed using Algorithm <ref> and each of the
π_t^l,N_l,M_l(φ) - π_t^l-1,N_l,M_l(φ) are produced by using Algorithm <ref>.
The cost of computing this estimator, per time step, is 

    𝒪(∑_l=0^LΔ_l^-1M_l(M_l+N_l)).

In the next section we shall establish how to choose L, (N_0,M_0),…,(N_L,M_L) so as to obtain a particular MSE and minimize the cost.










§ THEORETICAL ANALYSIS



 §.§ Notation

Define the sets ℕ={1,…} and ℕ_0={0,1,…}. Denote by ℬ_b(ℝ^d_1,ℝ^d_2) the set of all ℝ^d_2 valued bounded Borel measurable functions defined on ℝ^d_1, for k∈ℕ denote by 𝒞^k_b(ℝ^d_1,ℝ^d_2) the set of k-times continuously differentiable ℝ^d_2-valued functions whose domain is ℝ^d_1 and whose derivatives of order at most k with respect to any components are bounded. Define 𝒞^∞_b(ℝ^d_1,ℝ^d_2)=∩_k=1^∞𝒞^k_b(ℝ^d_1,ℝ^d_2). We use the shortcuts 𝒞_b^k(ℝ^d)=𝒞_b^k(ℝ^d,ℝ) and ℬ_b(ℝ^d)=ℬ_b(ℝ^d,ℝ).
For φ∈ℬ_b(ℝ^d) define the norm |φ|_0=x∈ℝ^dsup|φ(x)| and for φ∈𝒞_b^k(ℝ^d) define the norms 
    |φ|_k=∑_i=0^ksup_(j_1,…,j_i)∈{1,…,d_1}^isup_x∈ℝ^d|∂^i/∂ x_j_1…∂ x_j_iφ(x)|,

where x=(x_1,…,x_d_1). For a vector v∈ℝ^d denote by v the usual Euclidean norm.



 §.§ Assumptions


The functions a∈𝒞^2_b(ℝ^d×ℝ,ℝ^d)∩ℬ_b(ℝ^d×ℝ,ℝ^d), b∈𝒞^2_b(ℝ^d×ℝ,ℝ^d× d)∩ℬ_b(ℝ^d×ℝ,ℝ^d× d), and ξ_1,ξ_2∈𝒞^2_b(ℝ^d×ℝ^d)∩ℬ_b(ℝ^d×ℝ^d), and 
    inf_x∈ℝ^d+1inf_v∈ℝ^d\{0} v^⊤b(x)^⊤b(x)v/v^2 > 0.



The function G∈𝒞_b^2(ℝ^d)∩ℬ_b(ℝ^d) and satisfies inf_(x,y)∈ℝ^d× d_yG(x,y)>0.


	The sequence of positive integers N_0,…,N_L used to define the estimator π^L is strictly decreasing and satisfies N_l - N_l+1≥ L-l for 0≤ l≤ L-1.


For every function φ∈𝒞_b^∞(ℝ^d)∩ℬ_b(ℝ^d) the function x↦ P_μ_t-1,t(φ)=∫φ(z)P_μ_t-1,t(x,dz) belongs to 𝒞_b^∞(ℝ^d).


The functions a∈𝒞^∞_b(ℝ^d×ℝ,ℝ^d), b∈𝒞^∞_b(ℝ^d×ℝ,ℝ^d× d), and ξ_1,ξ_2∈𝒞^∞_b(ℝ^d×ℝ^d). Futhermore there exist functions ζ_1,ζ_2∈𝒞^∞_b(ℝ^d) such that ξ_i(x,y)=ζ_i(x-y) for all x,y∈ℝ^d and i∈{1,2}.




 §.§ Main Theorems


Assume [assump:A1](A1-2) and assump:W1(W1-2). Then for any φ∈𝒞^∞_b(ℝ^d)∩ℬ_b(ℝ^d) and t∈{1,…,T} there exists a constant C<+∞ such that for every (L,M,N)∈ℕ_0×ℕ×ℕ:
	
    𝔼[ (π_t(φ) - π_t^L,M,N(φ))^2 ] ≤ C(Δ_L^2+1/N + 1/M).




	Assume [assump:A1](A1-3) and assump:W1(W1-2). Then for any φ∈𝒞^∞_b(ℝ^d)∩ℬ_b(ℝ^d) and t∈{1,…,T} there exists C<+∞ such that for every L∈ℕ_0, and sequences of positive integers M_0,…,M_L and N_0,…,N_L:
	
    𝔼[ (π_t(φ)-π^L_t(φ))^2 ]≤ C(Δ_L^2 + ∑_l=0^LΔ_l^1/2/N_l + L∑_l=0^LΔ_l/M_l).


Theorems <ref> and <ref> are proven in the appendix. Assumptions [assump:W1](W1-2) are not needed to bound the variance of the estimates π_t^L,M,N(φ) and π_t^L(φ).




 §.§ Comments on Computational Cost


Consider the estimator π^L,M,M_t(φ), to achieve an error of order 𝒪(ϵ^2) for ϵ>0. Theorem <ref> provides the guideline of selecting L∝ |log(ϵ)| and M,N∝ϵ^-2, this choice leads to a cost of order 𝒪(M(M+N)Δ_L) = 𝒪(ϵ^-5). The multi-level estimator π_t^L(φ) can achieve the same order of error 𝒪(ϵ^2) at a lower cost. By utilizing the right hand side of the inequality in Theorem <ref> as a constraint for the error and optimizing the cost given this error constraint using the Lagrange multipliers method, we arrive to the following choices L∝ |log(ϵ)|, M_l∝ϵ^-2|log(ϵ)|Δ_l^5/6 and N_l∝ϵ^-2Δ_l^1/3, the resulting cost using these choices is 𝒪(ϵ^-4|log(ϵ)|^2).  In certain scenarios, we believe it is possible to reduce the cost to 𝒪(ϵ^-4). 
This conjecture should be
justified by results that show better bounds on the weak error of the output law of Algorithm <ref> <cit.>, this by itself is not enough, but one expects that there are scenarios
where the bound in Theorem <ref> can be improved to the reduced error 𝒪(Δ_L^2 + ∑_l=0^LΔ_l^1/2/N_l +∑_l=0^LΔ_l/M_l). This reduced error allows us to reduce the cost to 𝒪(ϵ^-4).
Assumption [assump:A3](A3) states that the number of particles should decrease relatively quickly. 
The increase rate stated in the assumption is very mild and not restrictive at all.  In fact N_0,…,N_L are expected to decrease exponentially which is the case for the optimal choice discussed above.





§ NUMERICAL RESULTS



 §.§ Models



  §.§.§ Kuramoto Model

We consider the Kuramoto model {X_t}_t≥0 defined by the Mckean-Valsov SDE

    dX_t = (θ + ∫sin(X_t - y)dμ_t(y))dt + σ dW_t,   X_0 = x_0∈ℝ,  t∈[0,T],

where μ_t is the law of X_t, θ be a random variable, σ>0, and T∈ℕ. The Kuramoto model is a well-known mean field game model and it is frequently used for numerical demonstrations. 
Note that we expect that we can obtain the improved cost of 𝒪(ϵ^-4) in this case, due to the fact that the diffusion coefficient is constant.



  §.§.§ Modified Kuramoto Model

In this demonstration consider the process X_t_t≥0 the Mckean-Vlasov SDE

    dX_t = (θ + ∫sin(X_t - y)dμ_t(y))dt + (σ/1+X_t^2)dW_t,   X_0 = x_0∈ℝ,  t∈[0,T],

where μ_t is the law of X_t, θ be a random variable, σ>0, and T∈ℕ. We call it the modified Kuramoto Model.  We consider this case as we expect that the non-constant diffusion coefficient
should demonstrate our theorerical results in the previous section, in that the computational cost should be
 𝒪(ϵ^-4|log(ϵ)|^2).



 §.§ Simulation Results

For both models we impose observations {Y_k}_t=1^T at unit times that follow Y_k|X_k=x_k ∼𝒩(x_k, τ^2) for k∈{1,…,T} and set x_0=1, θ=0, σ=0.2, τ=1, and T=50. We simulated a path from each model using Algorithm <ref> with level 9 and 5000 particles. We compute the estimator π^L,M,N_T(φ) and π^L_T(φ) for the function φ(x)=x where M,N are chosen as discussed in subsection <ref>. In order to calculate the mean square errors (MSE) of the two estimators we considered a proxy ground truth for the filter, this ground truth was obtained by computing π^L,M,N_T(φ) at level 7. The MSEs in our simulation were calculated using 128 independent replications. 

Figure <ref> shows log-log plots for the mean square errors against the theoretical costs for the two filter estimators applied to filtering the Kuramoto model and the modified Kuramoto.
In both cases we expect the multilevel methods to have rates of about -2 (Kuramot0) and -2.2,
reflecting the costs of 𝒪(ϵ^-4) and  𝒪(ϵ^-4|log(ϵ)|^2),
 and this is exactly what we observe.  For the case of the PF (single level) the rate should be about -2.5, which
gives the cost  𝒪(ϵ^-5).  The plots show the expected improvements of the multilevel method.



















  §.§.§ Acknowledgements


AJ was supported by start-up funding at CUHK-SZ.



§ PROOFS



 §.§ Structure


The following appendix contains several technical results that are used to prove Theorems <ref> and <ref>.  We begin by a short section on the additional notations and auxiliary processes that are to be used in the appendix. Section <ref> is devoted to studying the auxiliary processes {Z^l,M,i}_i=1^M, which are defined below. Sections <ref> and <ref> provide technical results leading to the proof of Lemma <ref>. Section <ref> introduces further several technical results and employs Lemma <ref> to ultimately prove Theorem <ref>. Section <ref> presents a technical discussion and uses the previously proven results in order to present a proof for Theorem <ref>. To understand the appendix, it should be read in chronological order.



 §.§ Notation

Let (Ω,ℱ,ℙ) be a measure space, W_t a Brownian motion on Ω, and X_t be the solution of SDE (<ref>). For random variables A_1,…,A_n and σ-subalgebras ℱ_1,…ℱ_n we define σ(A_1,…,A_n,ℱ_1,…,ℱ_n) to be the smallest σ-subalgebra 𝒢 of ℱ that contains the set ⋃_i=1^nℱ_i and that satisfies that the random variables A_1,…,A_n are measurable with respect to 𝒢. We consider step 2. of Algorithm <ref>: This step calls Algorithm <ref>, and step 2. of Algorithm <ref> is called T times. For (l,M,N)∈ℕ^3 given, we then use the associated notation:


    
  * {μ_kΔ_l^l,M}_k=0^k=TΔ_l^-1 and {μ_kΔ_l-1^l-1,M}_k=0^k=TΔ_l-1^-1 are the output measures of all the iterations of step 2. in <ref>.
    
  * {{Z_kΔ_l^l,M,i}_k=0^k=TΔ_l^-1}_i=1^M and {{Z_kΔ_l-1^l-1,M,i}_k=0^k=TΔ_l-1^-1}_i=1^M are the particles that generated the outputs of all the iterations of step 2. in <ref>. Let ℱ^l,M_kΔ_l = σ({{Z_jΔ_l^l,M,i}_j=0^j=k}_i=1^M) and ℱ^l-1,M_kΔ_l-1 =σ({{Z_jΔ_l-1^l-1,M,i}_j=0^j=k}_i=1^M) be the natural filterations of these particles respectively.
    
  * {B^i_t}_i=1^M are the Brownian motions generated in all the iterations of step 2. in Algorithm <ref> in order to generate the particles {{Z_kΔ_l^l,M,i}_k=1^k=TΔ_l^-1}_i=1^M and {{Z_kΔ_l-1^l-1,M,i}_k=1^k=TΔ_l-1^-1}_i=1^M. Denote Δ B^l,i_t = B^i_t - B^i_t-Δ_l for t∈[Δ_l,T] and i∈{1,…,M}. Similarly denote Δ W^l_t = W_t - W_t-Δ_l for t∈[Δ_l,T].
    
  * For s∈{l-1,l} define the discrete process {{Z_kΔ_s^s,i}_k=0^k=TΔ_l^-1}_i=1^M as the solution of the equations
        
        
        
        
    Z_(k+1)Δ_s^s,i  =  Z_kΔ_s^s,i + a(Z_kΔ_s^s,i,ξ_1(Z_kΔ_s^s,i,μ_kΔ_l^s)) + b(Z_kΔ_s^s,i,ξ_2(Z_kΔ_s^s,i,μ_kΔ_l^s))Δ B^s,i_(k+1)Δ_s,

        with μ^s_kΔ_s being the law of Z^s,1_kΔ_s and with the initial condition Z^s,i_0=x_0 for every (k,i)∈{0,…,TΔ_s^-1-1}×{1,…,M}.
    
    
  * For s∈{l-1,l} define the processes Z̃^s,M_t, and Z̃^s_t, as the solutions of 
            
    Z̃_(k+1)Δ_s^s  =  Z̃_kΔ_s^s + a(Z̃_kΔ_s^s,ξ_1(Z̃_kΔ_s^s,μ_kΔ_s^s)) +b(Z̃_kΔ_s^s,ξ_2(Z̃_kΔ_s^s,μ_kΔ_s^s))Δ W^s_(k+1)Δ_s,
    Z̃_(k+1)Δ_s^s,M =  Z̃_kΔ_s^s,M + a(Z̃_kΔ_s^s,M,ξ_1(Z̃_kΔ_s^s,M,μ_kΔ_s^s,M)) +b(Z̃_kΔ_s^s,M,ξ_2(Z̃_kΔ_s^s,M,μ_kΔ_s^s,M))Δ W^s_(k+1)Δ_s,

            with μ_kΔ_s^s being the law of Z_kΔ_s^s,1, μ_kΔ_s^s,M=1/M∑_i=1^Mδ_Z_kΔ_s^s,M,i and with the initial condition Z̃^s_0=Z̃^s,M_0=x_0 for every k∈{0,…,TΔ_s^-1-1}. Conditioning on the  σ-algebra ℱ^s,M_t denote by P_μ_t-1^s,M,t^s,M the transition kernel the process {Z̃^s,M_t}_t∈{0,1,…,T}.
    
    
  * For φ∈ℬ(ℝ^d) define
        
    γ_t(φ) := ∫_ℝ^dtφ(x_t){∏_p=1^t G(x_p,y_p)}∏_p=1^tP_μ_p-1,p(x_p-1,dx_p),

        
    γ_t^l(φ) := ∫_ℝ^dtφ(x_t){∏_p=1^t G(x_p,y_p)}∏_p=1^tP_μ_p-1^l,p^l(x_p-1,dx_p),

        and for s∈{l-1,l} define
        
    γ_t^l,M(φ) := ∫_ℝ^dtφ(x_t){∏_p=1   ^t G(x_p,y_p)}∏_p=1^tP_μ_p-1^l,M,p^l,M(x_p-1,dx_p).

    
  * Define π^l,M_t(φ) = γ^l,M_t(φ)/γ^l,M_t(1).


The Brownian motions {B_t^i}_i=1^M and W_t are independent. Notice that γ_t^l,M(φ) and γ_t^l-1,M(φ) are random variables measurable with respect to σ(ℱ^l-1,M_t,ℱ^l,M_t). The process Z̃_t^s is the Euler-Maruyama discretization for the process X_t. We introduce the process Z̃_t^s,M to approximate Z̃_t^s. The process Z̃_t^s,M utilizes the particles {Z_t^s,M,i}_i=1^M to create an empirical measure that approximates the measure μ^s_t independently of W_t. Our main objective is study the increments π_t^l,M,N(φ)-π_t^l-1,M,N(φ). To achieve this, we investigate the increments Z_t^l,M,i-Z_t^l-1,M,i and Z̃_t^l,M-Z̃_t^l-1,M. Finally we leverage the obtained results and the theory of multi-level particle filters to bound the 𝕃_2 error of the estimator π_t^L(φ).

Throughout the proofs below we will denote constants by C. For concise notation C may change values from step to another. Unless explicitly stated C is independent of l,M,N, and i.



 §.§ Bounds on Z^l,M,i_t increments



    Let Δ>0 and let {u_kΔ} be a sequence of non-negative real numbers indexed by the set Δℕ_0={kΔ:k∈ℕ_0}. Suppose that there exist constants C_1,C_2>0 for the which for every k∈ℕ the following inequality holds:
    
    u_kΔ≤ (1+C_1Δ)u_(k-1)Δ + C_2Δ.

    Then, the following inequality holds for every k∈ℕ_0
    
    u_kΔ≤ e^kΔ C_1(u_0 + kΔ C_2).



    Let k∈ℕ, there exists C>0 such that for all (A,B)∈(ℝ^d)^2 The following inequality holds
    
    A+B^2k≤A^2k + 2kA^2k-2A^TB + CB^2(A^2k-2+B^2k-2).




Assume [assump:A1](A1). Let k∈ℕ. There exists a constant C<+∞ such that for every (l,M)∈ℕ_0×ℕ, t∈{0,Δ_l,…,T}, and i∈{1,…,M}:
    
    max{𝔼[Z^l,M,i_t^k], 𝔼[Z^l,i_t^k]}≤ C,



    It is enough to prove the inequality for even positive integers k, the general case follows from the inequality 𝔼[A^k]≤𝔼[A^2k]^1/2 for any ℝ^d valued random variable A. Using the definition (<ref>) of Z^l,i_t, applying Lemma <ref> twice, Assumption [assump:A1](A1) that the functions a,b are bounded, Δ_l≤ 1, 𝔼[Δ B^l,i_t-Δ_l]=0, and 𝔼[B^l,i_t-Δ_l^k]≤ CΔ_l^k/2 we have
    
    𝔼[Z^l,i_t^k] =    𝔼[ Z^l,i_t-Δ_l + a(Z^l,i_t-Δ_l,ξ̅_2(Z^l,i_t-Δ_l,μ^l_t-Δ_l))Δ + b(Z^l,i_t-Δ_l,ξ̅_2(Z^l,i_t-Δ_l,μ^l_t-Δ_l))Δ B^l,i_t-Δ_l^k]
    ≤    𝔼[Z^l,i_t-Δ_l^k] + CΔ_lmax_r∈{1,2,3,4}𝔼[Z^l,i_t-Δ_l^k-r] + CΔ_l^k/2
    ≤     (1+CΔ_l)𝔼[Z^l,i_t-Δ_l^k] + CΔ_l.

    The last line follows from the inequality x^s≤ 1+x^k for 0≤ s≤ k and x≥0.
    Applying Lemma <ref> to the sequence {[Z^l,i_sΔ_l^k]}_s=0^T/Δ_l yields the inequality 𝔼[Z^l,i_t-Δ_l^k]≤ C, clearly following similar steps yields 𝔼[Z^l,M,i_t-Δ_l^k]≤ C. 



    Assume [assump:A1](A1). Let k∈ℕ. There exists a constant C<+∞ such that for every (l,M)∈ℕ×ℕ, t∈{0,Δ_l,…,T}, and i∈{1,…,M}:
    
    𝔼[ Z^l,M,i_t-Z^l,i_t^k ]≤ C1/M^k/2,

    
    𝔼[ Z^l,M,i_t-Z^l-1,M,i_t^k ]≤ CΔ_l^k/2,



    Let t∈{Δ,2Δ,…,T}, we have.
	
    Z^l,M,i_t - Z^l,i_t =     Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l
    +     (a(Z^l,M,i_t-Δ_l,ξ̅_1(Z^l,M,i_t-Δ_l,μ^l,M_t-Δ_l))-a(Z^l,i_t-Δ_l,ξ̅_1(Z^l,i_t-Δ_l,μ^l_t-Δ_l)))Δ_L 
     +     (b(Z^l,M,i_t-Δ_l,ξ̅_2(Z^l,M,i_t-Δ_l,μ^l,M_t-Δ_l))-b(Z^l,i_t-Δ_l,ξ̅_2(Z^l,i_t-Δ_l,μ^l_t-Δ_l)))Δ B^l,i_t.

    Using the C^p inequality and that the functions a and ξ_1 are bounded and Lipchitz (Assumption [assump:A1](A1)) there exists a constant C>0 that satisfies

    𝔼[a(Z^l,M,i_t-Δ_l,ξ̅_1(Z^l,M,i_t-Δ_l,μ^l,M,i_t-Δ_l))-a(Z^l,i_t-Δ_l,ξ̅_1(Z^l,i_t-Δ_l,μ^l,i_t-Δ_l))^k]
    ≤     C(𝔼[Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l^k]  + 𝔼[|1/M∑_j=1^Mξ_1(Z^l,M,i_t-Δ_lZ^l,M,j_t-Δ_l)-1/M∑_j=1^Mξ_1(Z^l,i_t-Δ_l,Z^l,j_t-Δ_l)|^k]
    
            +    1/M^k𝔼[|∑_j=1
    j≠i^M(ξ_1(Z^l,i_t-Δ_l,Z^l,j_t-Δ_l)-∫ξ_1(Z^l,i_t-Δ_l,y)dμ^l_t-Δ_l(y))|^k]  ) + C1/M^k.

    The second term after the inequality is bounded by 1/M∑_j=1^M𝔼[|ξ_1(Z^l,M,i_t-Δ_lZ^l,M,j_t-Δ_l)-ξ_1(Z^l,i_t-Δ_l,Z^l,j_t-Δ_l)|^k]≤ Csup_i𝔼[Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l^k] by the C^p inequality. 
    For the third term we use the idependence of Z_t^l,i and Z_t^l,j for i≠j, Fubini Theorem, and Marcinkiewicz-Zygmund inequality: 
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    
    𝔼[|∑_j=1
    j≠i^M(ξ_1(Z^l,i_t-Δ_l,Z^l,j_t-Δ_l)-∫ξ_1(Z^l,i_t-Δ_l,y)dμ^l_t-Δ_l(y))|^k ]
    
                 =    ∫𝔼[|∑_j=1
    j≠i^M(ξ_1(z,Z^l,j_t-Δ_l)-∫ξ_1(z,y)dμ^l_t-Δ_l(y)d)|^k ]μ^l_t-Δ_l(z)
    ≤     C∫𝔼[ (∑_j=1
    j≠i^M(ξ_1(z,Z^l,j_t-Δ_l)-∫ξ_1(z,y)dμ^l_t-Δ_l(y))^2)^k/2]μ^l_t-Δ_l(z)
    ≤     CM^k/2,

    where the last line follows from the boundedness of ξ_1. Thus we have shown that
    
    𝔼[a(Z^l,M,i_t-Δ_l,ξ̅_1(Z^l,M,i_t-Δ_l,μ^l,M,i_t-Δ_l))-a(Z^l,i_t-Δ_l,ξ̅_1(Z^l,i_t-Δ_l,μ^l,i_t-Δ_l))^k]
                ≤ Csup_i𝔼[Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l^k] + C1/M^k/2.

    A similar bound holds for the function b. Analogously to the proof of Lemma <ref>, we apply Lemma <ref> twice on Z^l,M,i_t-Z^l,i_t following equation (<ref>) and using the bound (<ref>) (for both a and b) yields
    
    sup_i𝔼[Z^l,M,i_t - Z^l,i_t^k]≤ (1+CΔ_l)sup_i𝔼[Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l^k] + CΔ_l1/M^k/2.

    Applying Lemma <ref> to the sequence {sup_i𝔼[Z^l,M,i_sΔ_l - Z^l,i_sΔ_l^k]}_s=0^T/Δ_l and using Z^l,i_0=Z^l,M,i_0=x_0 proves the first inequality. For the second inequality we write
    
    Z^l,M,i_t  - Z^l-1,M,i_t
    
                 =     Z^l,M,i_t-Δ_l-1  - Z^l-1,M,i_t-Δ_l-1
    
                 +    (a(Z^l,M,i_t-Δ_l-1,ξ̅_1(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)) - a(Z^l-1,M,i_t-Δ_l-1,ξ̅_1(Z^l-1,M,i_t-Δ_l-1, μ^l-1,M_t-Δ_l-1)))Δ_l-1
    
                 +    (b(Z^l,M,i_t-Δ_l-1,ξ̅_2(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)) - b(Z^l-1,M,i_t-Δ_l-1,ξ̅_2(Z^l-1,M,i_t-Δ_l-1, μ^l-1,M_t-Δ_l-1)))Δ B^l-1,i_t
    
                 +    (a(Z^l,M,i_t-Δ_l,ξ̅_1(Z^l,M,i_t-Δ_l, μ^l,M_t-Δ_l)) - a(Z^l,M,i_t-Δ_l-1,ξ̅_1(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)))Δ_l
    
                 +    (b(Z^l,M,i_t-Δ_l,ξ̅_2(Z^l,M,i_t-Δ_l, μ^l,M_t-Δ_l)) - b(Z^l,M,i_t-Δ_l-1,ξ̅_2(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)))Δ B^l,i_t.

    Bounding each term following the same technique used to prove the first inequality we can show the following the inequality
        
    sup_i𝔼[Z^l,M,i_t  - Z^l-1,M,i_t^k]
    ≤     (1+CΔ_L)sup_i𝔼[Z^l,M,i_t-Δ_l-1  - Z^l-1,M,i_t-Δ_l-1^k]
                + CΔ_l-1(sup_i𝔼[Z^l,M,i_t-Δ_l  - Z^l,M,i_t-Δ_l-1^k] + Δ_l-1^k/2).

    From the definition of Z^l,M,i_t and the boundedness of a and b we have sup_i𝔼[Z^l,M,i_t-Δ_l  - Z^l,M,i_t-Δ_l-1^k]≤ CΔ_l^k/2. The proof then follows similarly to the proof of the first inequality by applying Lemma <ref> to the sequence {sup_i𝔼[Z^l,M,i_sΔ_l-1  - Z^l-1,M,i_sΔ_l-1^k]}_s=0^T/Δ_l-1.






    Let f∈𝒞_b^2(ℝ^p,ℝ^q). There exists a constant C<+∞ such that the following inequality holds for every (x,y,z,w)∈(ℝ^p)^4:
    
    f(x)-f(y)-f(z)+f(w)≤ Cx-y-z+w + Cz-w(z-x+w-y).





    Assume [assump:A1](A1). Let k∈ℕ. There exists a constant C<+∞ such that for every (l,M)∈ℕ×ℕ, t∈{0,Δ_l,…,T}, and i∈{1,…,M}:
    
    𝔼[ Z^l,M,i_t-Z^l-1,M,i_t - Z^l,i_t+Z^l-1,i_t^k ]≤ CΔ_l^k/2/M^k/2.



Assume k is even. The inequality is proven following the same strategy as the proof of Lemma <ref>, we write following (<ref>) and (<ref>):
    
    Z^l,M,i_t  - Z^l-1,M,i_t -  Z^l,i_t  + Z^l-1,i_t
    
        		=     Z^l,M,i_t-Δ_l-1  - Z^l-1,M,i_t-Δ_l-1 - Z^l,i_t-Δ_l-1  + Z^l-1,i_t-Δ_l-1
    
        		+    (a(Z^l,M,i_t-Δ_l-1,ξ̅_1(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)) - a(Z^l-1,M,i_t-Δ_l-1,ξ̅_1(Z^l-1,M,i_t-Δ_l-1, μ^l-1,M_t-Δ_l-1))
    
        		-     a(Z^l,i_t-Δ_l-1,ξ̅_1(Z^l,i_t-Δ_l-1, μ^l_t-Δ_l-1)) + a(Z^l-1,i_t-Δ_l-1,ξ̅_1(Z^l-1,i_t-Δ_l-1, μ^l-1_t-Δ_l-1)))Δ_l-1
    
        		+    ⋯

where we wrote ⋯ for brevity to indicate that the rest of terms are similar to first term in (<ref>) and follows the same general form as (<ref>) and (<ref>). We will show how to bound the first term in (<ref>) since the rest of the terms, which were indicated by ⋯, are treated similarly. define the vector R^l,M = (Z^l,M,i_t-Δ_l-1,ξ̅_1(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)) and write a(R^l,M) for a(Z^l,M,i_t-Δ_l-1,ξ̅_1(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)), similarly define R^l-1,M, R^l, and R^l-1. Applying Lemma <ref> to the function a we have
    
    a(R^l,M) - a(R^l-1,M) - a(R^l) + a(R^l-1)≤    CR^l,M - R^l-1,M(R^l,M - R^l+R^l-1,M - R^l-1) 
    
                +     CR^l,M - R^l-1,M - R^l + R^l-1,

    Using the assumption that ξ_1 is Liptchitz, Lemma <ref>, and Marcinkiewicz-Zygmund in a manner similar to (<ref>) we have 𝔼[R^l,M-R^l-1,M^2k]≤  CΔ_l^k and 𝔼[R^l,M-R^l^2k]≤ C1/M^k,
    thus using Cauchy-Schwarz and the C^p inequality we have
    
    𝔼[R^l,M - R^l-1,M^k(R^l,M - R^l+R^l-1,M - R^l-1)^k]
    ≤    C(𝔼[R^l,M - R^l,M^2k](𝔼[R^l,M - R^l^2k]+𝔼[R^l-1,M - R^l-1^2k]))^1/2
    ≤     CΔ_l^k/2/M^k/2.

    Similarly to (<ref>) the term R^l,M - R^l-1,M - R^l + R^l-1 in (<ref>) is bounded by
    
    𝔼[R^l,M - R^l-1,M - R^l + R^l-1^k]
    ≤     C𝔼[Z_t-Δ_l-1^l,M,i - Z_t-Δ_l-1^l-1,M,i - Z_t-Δ_l-1^l,i + Z_t-Δ_l-1^l-1,i^k]
    
                +     C𝔼[|1/M∑_j=1^M(ξ_1(Z^l,M,i_t-Δ_l-1, Z^l,M,j_t-Δ_l-1) - ξ_1(Z^l-1,M,i_t-Δ_l-1, Z^l-1,M,j_t-Δ_l-1)
                -ξ_1(Z^l,i_t-Δ_l-1, Z^l,j_t-Δ_l-1) 
    
                +    ξ_1(Z^l,i_t-Δ_l-1, Z^l,j_t-Δ_l-1))|^k]
                + C𝔼[|1/M∑_j≠i(ξ_1(Z^l,i_t-Δ_l-1, Z^l,j_t-Δ_l-1) - ξ_1(Z^l-1,i_t-Δ_l-1, Z^l-1,j_t-Δ_l-1)
    
                -    ∫ξ_1(Z^l,i_t-Δ_l-1, y)dμ^l_t-Δ_l-1(y) + ∫ξ_1(Z^l-1,i_t-Δ_l-1, y)dμ^l-1_t-Δ_l-1(y) )|^k]
                + C1/M^k𝔼[|ξ_1(Z^l,i_t-Δ_l-1, Z^l,i_t-Δ_l-1)
     -    ξ_1(Z^l-1,i_t-Δ_l-1, Z^l-1,i_t-Δ_l-1)
                -∫ξ_1(Z^l,i_t-Δ_l-1, y)dμ^l_t-Δ_l-1(y) + ∫ξ_1(Z^l-1,i_t-Δ_l-1, y)dμ^l-1_t-Δ_l-1(y)|^k ].

    Using Lemmata <ref> and then <ref> and with Cauchy-Schwarz similarly to <ref>, the second expectation after the inequality in (<ref>) can be bounded from above by Csup_i 𝔼[Z_t-Δ_l-1^l,M,i - Z_t-Δ_l-1^l-1,M,i - Z_t-Δ_l-1^l,i + Z_t-Δ_l-1^l-1,i^k] + CΔ_L^k/2/M^k/2. For the third expectations we use Fubini and Marcinkiewicz-Zygmund inequality in a manner similar to (<ref>), ξ_1 is Lipchitz, and Lemma <ref>:
    
    𝔼[|∑_j≠i(ξ_1(Z^l,i_t-Δ_l-1, Z^l,j_t-Δ_l-1) - ξ_1(Z^l-1,i_t-Δ_l-1, Z^l-1,j_t-Δ_l-1)
                -∫ξ_1(Z^l,i_t-Δ_l-1, y)dμ^l_t-Δ_l-1(y)
     +    ∫ξ_1(Z^l-1,i_t-Δ_l-1, y)dμ^l-1_t-Δ_l-1(y) )|^k]
    ≤     C𝔼[|∑_j≠i(ξ_1(Z^l,i_t-Δ_l-1, Z^l,j_t-Δ_l-1) - ξ_1(Z^l-1,i_t-Δ_l-1, Z^l-1,j_t-Δ_l-1)
                -𝔼[ξ_1(Z^l,i_t-Δ_l-1, Z^l,j_t-Δ_l-1)]
     +    𝔼[ξ_1(Z^l-1,i_t-Δ_l-1, Z^l-1,j_t-Δ_l-1)] )^2|^k/2]
    ≤     CM^k/2sup_i 𝔼[Z^l,i_t-Δ_l-1-Z^l,i_t-Δ_l-1^k]
    ≤     CM^k/2Δ_l^k/2.

    Using ξ_1 is Lipchitz and Lemma <ref> the last expectation in (<ref>) is bounded by CΔ_l^k/2/M^k. Therefore we have shown
    
    𝔼[a(Z^l,M,i_t-Δ_l-1,ξ̅_1(Z^l,M,i_t-Δ_l-1, μ^l,M_t-Δ_l-1)) - a(Z^l-1,M,i_t-Δ_l-1,ξ̅_1(Z^l-1,M,i_t-Δ_l-1, μ^l-1,M_t-Δ_l-1))
    
        		-     a(Z^l,i_t-Δ_l-1,ξ̅_1(Z^l,i_t-Δ_l-1, μ^l_t-Δ_l-1)) + a(Z^l-1,i_t-Δ_l-1,ξ̅_1(Z^l-1,i_t-Δ_l-1, μ^l-1_t-Δ_l-1))^k]
    ≤     Csup_i𝔼[Z^l,M,i_t-Δ_l-1  - Z^l-1,M,i_t-Δ_l-1 -  Z^l_t-Δ_l-1  + Z^l-1_t-Δ_l-1^k]
                + CΔ_l^k/2/M^k/2.

    A similar bound holds for the rest of terms in (<ref>). Applying Lemma <ref> several times we have
    
    sup_i𝔼[Z^l,M,i_t  - Z^l-1,M,i_t -  Z^l_t  + Z^l-1_t^k]
    ≤     (1+CΔ_l)sup_i𝔼[Z^l,M,i_t-Δ_l-1  - Z^l-1,M,i_t-Δ_l-1 -  Z^l_t-Δ_l-1  + Z^l-1_t-Δ_l-1^k]
                + CΔ_l^k/2/M^k/2.

    Finally applying Lemma <ref> to the sequence {sup_i𝔼[Z^l,M,i_sΔ_l-1  - Z^l-1,M,i_sΔ_l-1 -  Z^l,i_sΔ_l-1  + Z^l-1,i_sΔ_l-1^k]}_s=0^TΔ_l-1^-1 finishes the proof.    



 §.§ Bounds on Z̃^l,M_t and π^l,M,N_t increments



Assume [assump:A1](A1). Let k∈ℕ. There exists a constant C<+∞ such that for every (l,M)∈ℕ×ℕ and t∈{0,Δ_l,…,T} the following inequalities hold
    
    max{𝔼[Z̃^l,M_t^k], 𝔼[Z̃^l_t^k]}≤ C,

	
    𝔼[Z̃^l,M_t - Z̃^l_t^k]≤ C1/M^k/2,
 
    
    𝔼[Z̃^l,M_t - Z̃^l-1,M_t^k]≤ CΔ_L^k/2 ,

	
    𝔼[Z̃^l,M_t - Z̃^l-1,M_t - Z̃^l_t + Z̃^l-1_t^k]≤ CΔ_l^k/2/M^k/2.



    We assume k is even. We will prove the second inequality only since the same technique we will use to prove it and the techniques used for proving Lemmata <ref> and <ref> can directly be applied to prove the other inequalities in this lemma statement. We write
	
    Z̃^l,M_t - Z̃^l_t =    Z̃^l,M_t-Δ_l - Z̃^l_t-Δ_l
    +     (a(Z̃^l,M_t-Δ_l,ξ̅_1(Z̃^l,M_t-Δ_l,μ^l,M_t-Δ_l))-a(Z̃^l_t-Δ_l,ξ̅_1(Z̃^l_t-Δ_l,μ^l_t-Δ_l)))Δ_l 
     +     (b(Z̃^l,M_t-Δ_l,ξ̅_2(Z̃^l,M_t-Δ_l,μ^l,M_t-Δ_l))-b(Z̃^l_t-Δ_l,ξ̅_2(Z̃^l_t-Δ_l,μ^l_t-Δ_l)))Δ W^l_t.

    The functions a and ξ_1 are Lipschitz from Assumption [assump:A1](A1), hence for a constant C we have

    𝔼[a(Z̃^l,M_t-Δ_l,ξ̅_1(Z̃^l,M_t-Δ_l,μ^l,M_t-Δ_l))-a(Z̃^l_t-Δ_l,ξ̅_1(Z̃^l_t-Δ_l,μ^l_t-Δ_l))^k]
    ≤     C(𝔼[Z̃^l,M_t-Δ_l - Z̃^l_t-Δ_l^k] + 𝔼[|1/M∑_i=1^Mξ_1(Z̃^l,M_t-Δ_l,Z^l,M,i_t-Δ_l)-1/M∑_i=1^Mξ_1(Z̃^l_t-Δ_l,Z^l,i_t-Δ_l)|^k]
     
            +    𝔼[|1/M∑_i=1^Mξ_1(Z̃^l_t-Δ_l,Z^l,i_t-Δ_l)-∫ξ_1(Z̃^l_t-Δ_l,y)dμ^l_t-Δ_l(y)|^k].

The second expectation after the inequality in (<ref>) is bounded by C(sup_i𝔼[Z̃^l,M_t-Δ_l - Z̃^l_t-Δ_l^k] + 𝔼[Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l^k])
and from Lemma <ref> we have sup_i 𝔼[Z^l,M,i_t-Δ_l - Z^l,i_t-Δ_l]≤ C/M^k/2. The third expectation after the inequality in (<ref>) is bounded above by C/M^k/2 using Fubini and Marcinkiewicz-Zygmund inequality in a manner similar to (<ref>). Therefore we have

    𝔼[a(Z̃^l,M_t-Δ_l,ξ̅_1(Z̃^l,M_t-Δ_l,μ^l,M_t-Δ_l))-a(Z̃^l_t-Δ_l,ξ̅_1(Z̃^l_t-Δ_l,μ^l_t-Δ_l))^k]
    		≤  C(𝔼[Z̃^l,M_t-Δ_l - Z̃^l_t-Δ_l^k] + 1/M^k/2).

The rest of the proof is identical to the proof of Lemma <ref>.



	Assume [assump:A1](A1-2). Let φ∈𝒞^2_b(ℝ^d)∩ℬ_b(ℝ^d), t∈{0,1,2,…,T}, and k∈ℕ. There exists C<+∞ such that for every (l,M)∈ℕ×ℕ:
	
	
	
	
	
	
 	
    𝔼[|(γ_t^l,M - γ_t^l-1,M)(φ)|^k] ≤ C Δ_l^k/2|φ|_1^k,    𝔼[|(γ_t^l,M - γ_t^l)(φ)|^k] ≤ C 1/M^k/2|φ|_1^k,


	
    𝔼[|(γ_t^l,M - γ_t^l-1,M - γ_t^l + γ_t^l-1)(φ)|^k] ≤ C Δ_l^k/2/M^k/2|φ|_2^k.



    Let t∈{0,1,…,T}. Define the function Φ:(ℝ^d)^t→ℝ by Φ(x_1,…,x_t) = φ(x_t)∏_i=1^t G(x_i,y_i). From Assumption [assump:A2](A2) the function Φ is Lipchitz and |Φ|_1≤ C|φ|_1. Using Jensen's inequality, the mean value theorem, and Lemma <ref> we have
    
    𝔼[|(γ_t^l,M - γ^l-1,M_t)(φ)|^k] =    𝔼[|𝔼[ Φ(Z̃_1^l,M,…,Z̃_t^l,M) - Φ(Z̃_1^l-1,M,…,Z̃_t^l-1,M)| ℱ^l,M_t,ℱ^l-1,M_t ]|^k ]
    ≤    𝔼[ |Φ(Z̃_1^l,M,…,Z̃_t^l,M) - Φ(Z̃_1^l-1,M,…,Z̃_t^l-1,M)|^k ]
    ≤     C𝔼[ |φ|_1^k ∑_i=1^tZ̃_i^l,M - Z̃_i^l-1,M^k ]
    ≤     C|φ|_1^k Δ_l^k/2.

    The second inequality in the lemma statement follows analogously. The proof of the third inequality is essentially the same but we use Lemma <ref> and <ref>.


	Assume [assump:A1](A1-2). Let φ∈𝒞^2_b(ℝ^d)∩ℬ_b(ℝ^d), t∈{0,1,2,…,T}, and k∈ℕ. There exists C<+∞ such that for every (l,M)∈ℕ×ℕ:
	
	
	
	
	
	
 	
    𝔼[|(π_t^l,M - π_t^l-1,M)(φ)|^k] ≤ C Δ_l^k/2|φ|_1^k,    𝔼[|(π_t^l,M - π_t^l)(φ)|^k] ≤ C 1/M^k/2|φ|_1^k,


	
    𝔼[|(π_t^l,M - π_t^l-1,M - π_t^l + π_t^l-1)(φ)|^k] ≤ C Δ_l^k/2/M^k/2|φ|_2^k.



The first two inequalities follow directly form Lemma <ref>. The third inequality follows directly as well from Lemma <ref> and Cauchy Schwarz inequality and using the algebraic identity




    a/A - b/B - c/C + d/D =    1/A(a-b-c+d) -b/AB(A-B-C+D) - 1/AC(A-C)(c-d)
    
        -    1/AB(C-D)(b-d)+d/CBD(B-D)(C-D) + d/ACB(A-C)(C-D).






 §.§ Proof of Theorem <ref>


We begin this section by bounding the bias of the estimator π^L(φ). We use Proposition 3.1 from <cit.> which we restate as Lemma <ref> below.

    Consider the process Z_t' =Z_t-x_0 and define the functions a'(u,v)=a(u+x_0,v), b'(u,v)=b(u+x_0,v), ξ_j'(u,v)=ξ_j(u+x_0,v+x_0) for j∈{1,2}. Z_t' satisfies SDE (<ref>) but with a,b,ξ_1,ξ_2,x_0 replaced with a',b',ξ_1',ξ_2',0. C in the Proposition 3.1 in <cit.> can be bounded by a constant that depends on the bounds of a',b',ξ_1',ξ_2' and their derivatives which are finite by Assumption [assump:A1]A1, this allows us to remove the dependence of x_0 from the constant C in Proposition 3.1 of <cit.>.


	Assume [assump:W2](W2). Let t∈{0,1,…,T}. For every φ∈𝒞^∞_b(ℝ^d)∩ℬ_b(ℝ^d) there exists C<+∞ such that for every l∈ℕ_0:
     
    |P_μ_t-1,t(φ) - P_μ_t-1^l,t(φ)|_0≤ CΔ_l.



Assume [assump:A2](A2). Let t∈{0,…,T} There exists a constant C<+∞ such that for every l∈ℕ_0 and  φ∈ℬ_b(ℝ^d):
    
    max(|γ_t(φ)|,|γ_t^l(φ)|) ≤ C|φ|_0,



    Follows directly from Assumption [assump:A2](A2) and the definitions of γ_t(φ) and γ_t^l(φ).


Assume [assump:A1](A1-2) and [assump:W1](W1-2). Let t∈{0,…,T}. For every φ∈𝒞^∞_b(ℝ^d)∩ℬ_b(ℝ^d) there exists a constant C<+∞ such that for every (l,M)∈ℕ_0×ℕ
    
    |γ_t(φ) - γ_t^l(φ)|_0 ≤ CΔ_l,

    
    |π_t(φ) - π_t^l(φ)|_0 ≤ CΔ_l.



    For concise notation we define the functions G_t(x)=G(x,y_t) for every t∈{0,…,T} and x∈ℝ^d. We will prove the inequalities using induction.
    The first inequality is satisfied trivially at t=0 since γ_0=γ_0^L=δ_x_0. Assume the inequality holds for t-1 where t≥ 1, we write
    
    γ_t(φ) - γ_t^l(φ)
            =γ_t-1^l(P_μ_t-1,t(G_tφ) - P^l_μ^l_t-1,t(G_tφ)) + (γ_t-1-γ_t-1^l)(P_μ_t-1,t(G_tφ)).

    Lemma <ref> and Lemma <ref> yields the bound for the first term. Assumption [assump:W1](W1) allows us to apply the induction hypothesis to bound the second term. 
    
    
    
    The second inequality follows from the first one and the equality
    
    π_t(φ) - π_t^l(φ)=
                1/γ_t(1)γ_t^l(1)(γ_t^l(1)-γ_t(1))γ_t^l(φ) - 1/γ_t(1)(γ^l_t(φ) - γ_t(φ)).

    Notice that the quantities γ_t(1) and γ^l_t(1) are bounded below away from 0 because of Assumption [assump:A2](A2).

The following lemma is a consequence of the standard theory of 𝕃_p-bounds for particle filters (Theorem 2.9 in <cit.>). 





	Assume [assump:A1](A1-2) and let t∈{1,…,T}. There exists C<+∞ such that for every φ∈ℬ_p(ℝ^d) and (l,M,N)∈ℕ_0×ℕ×ℕ:
	
    𝔼[ |π^l,M_t(φ)- π^l,M,N_t(φ) |^2 ] ≤ C1/N|φ|_0^2.



    Because G is bounded we can apply Theorem 2.9 in <cit.> to obtain the bound 𝔼[ (π^l,M_t(φ)- π^l,M,N_t(φ))^2 |ℱ^l,M_t]≤ C|φ|_0/N where C is ℱ^l,M_t. However, following the formula for C in Theorem 2.9 in <cit.> and Assumption [assump:A2](A2) allow us to bound C from above by a deterministic constant. This justifies the deterministic C in Lemma <ref>.



Using the C^p inequality
    
    𝔼[(π_t(φ) - π_t^L,M,N(φ))^2]
    ≤     C( (π_t(φ) - π_t^L(φ))^2+ 𝔼[(π_t^L(φ) - π_t^L,M(φ))^2] + 𝔼[(π_t^L,M(φ) - π_t^L,M,N(φ))^2] ).

    Lemmata <ref>, <ref>, and <ref> finish the proof.



 §.§ Proof of Theorem <ref>


We condition on σ(ℱ^l-1,M_t,ℱ^l,M_t) and consider Theorem C.1 and Lemma C.3 in <cit.> under this conditioning. Following the same reasoning as in Remark <ref> the constant C(n,φ) in <cit.> can be assumed to be deterministic in our setting because the functions a,b,ξ_1,ξ_2,G are bounded by Assumptions [assump:A1](A1-2). Thus we have that for s∈{1,2}:

    |𝔼[ ( π^l,M,N_t(φ) - π^l-1,M,N_t(φ) - π^l,M_t(φ) + π^l-1,M_t(φ))^s]| ≤ C𝔼[B(n)]^s/2/N,

where we used the inequality 𝔼[√(B(n))]≤𝔼[B(n)]^1/2 and B(n) is defined in equation (C.1) in <cit.>. Lemma <ref> verifies an in expectation (over all randomness) version of Assumption D.1 in <cit.> from which a version of Theorem D.5 in <cit.> follows (with expectation taken over all randomness). We apply the C^p inequality on B(n), the first sum corresponding to maximum coupling is bounded by our version of Theorem D.5 in <cit.> and for the other terms Lemmata <ref> and <ref> suffice. Therefore we have 𝔼[B(n)]≤Δ^1/2. We state the result of this discussion as the following Lemma.


	Assume [assump:A1](A1-2). Let φ∈𝒞_b^1(ℝ^d)∩ℬ_b(ℝ^d) and t∈{1,…,T}. There exists C<+∞ such that for every (l,M,N)∈ℕ×ℕ×ℕ and s∈{1,2}:
	
    |𝔼[ ( π^l,M,N_t(φ) - π^l-1,M,N_t(φ) - π^l,M_t(φ) + π^l-1,M_t(φ))^s]| ≤ CΔ_l^s/4/N




    
    𝔼[ (π_t(φ)-π^L_t(φ))^2 ]
    
                =    𝔼[ ((π_t(φ) - π_t^L(φ)) + (π_t^0(φ) - π_t^0,M_0(φ)) + ( π_t^0,M_0(φ) - π_t^0,M_0,N_0(φ)) 
    
                +    ∑_l=1^L(π^l_t - π^l-1_t-π^l,M_l_t +π^l-1,M_l_t)(φ)  + ∑_l=1^L(π^l,M_l_t - π^l-1,M_l_t-π^l,M_l,N_l_t +π^l-1,M_l,N_l_t)(φ))^2 ]
    ≤    C(π_t(φ)-π_t^L(φ))^2 + C𝔼[ (π_t^0(φ) - π_t^0,M_0(φ)) ^2 ] + C𝔼[  ( ∑_l=1^L(π^l_t - π^l-1_t-π^l,M_l_t +π^l-1,M_l_t)(φ))^2]
    
                +     C𝔼[ ( π_t^0,M_0(φ) - π_t^0,M_0,N_0(φ))^2] + C𝔼[(∑_l=1^L(π^l,M_l_t - π^l-1,M_l_t-π^l,M_l,N_l_t +π^l-1,M_l,N_l_t)(φ))^2].

    The first term is bounded using CΔ_L^2 by Lemma <ref>, the bounds for the second and third terms follow from C^p inequality and Lemma <ref>. The fourth is the error of a standard particle filter stated in Lemma <ref>. For the fifth term we use Lemma <ref> as follows
    
    
    𝔼[(  ∑_l=1^L(π_t^l,M_l - π_t^l-1,M_l-π_t^l,M_l,N_l +π_t^l-1,M_l,N_l)(φ)   )^2]
    
                =    ∑_l=1^L𝔼[(  (π_t^l,M_l - π_t^l,M_l-π_t^l-1,M_l,N_l +π_t^l-1,M_l,N_l)(φ)    )^2]
                +∑_p≠q
        𝔼[ (π_t^p,M_p - π_t^p-1,M_p-π_t^p,M_p,N_p +π_t^p-1,M_p,N_p)(φ)  ]𝔼[ (π_t^q,M_q - π_t^q-1,M_q-π_t^q,M_q,N_q +π_t^q-1,M_q,N_q)(φ) ]
    ≤     C∑_l=1^LΔ_l^1/2/N_l + C∑_p≠qΔ_p^1/4Δ_q^1/4/N_pN_q
    ≤     C∑_l=1^LΔ_l^1/2/N_l,

	where for the last line first notice that Assumtion [assump:A3](A3) implies that N_l≥ (L-l)^2/2 for 0≤ l≤ L-1 then using Cauchy Schwarz inequality
	
    C∑_p≠qΔ_p^1/4Δ_q^1/4/N_pN_q≤ C(∑_l=1^LΔ_l^1/4/N_l)^2≤ C(∑_l=1^L1/N_l) (∑_l=1^LΔ_l^1/2/N_l)< C(1+1/2∑_l=1^∞1/l^2) (∑_l=1^LΔ_l^1/2/N_l).



99

bain
Bain, A., & Crisan, D. (2009). Fundamentals of Stochastic Filtering. Springer: New York.

barbu
Barbu,  V. & Roeckner, M. (2020). From nonlinear Fokker-Planck equations to solutions of distribution
dependent SDE.  Ann. Probab., 48, 1902–1920.

ml_is_mv
Ben Rached, N. , Haji-Ali, A.-L., Shyam, M., & Tempone, R. (2022).
Multilevel Importance Sampling for McKean-Vlasov Stochastic Differential Equation. arXiv preprint.

mi_is_mv
Ben Rached, N. , Haji-Ali, A.-L., Shyam, M., & Tempone, R. (2023).
Multi-Index Importance Sampling for McKean-Vlasov Stochastic Differential Equation.  arXiv preprint.

bush
Bush, N. , Hambly,  B. M.,  Haworth, H.,  Jin, L. & Reisinger, C.  (2011). Stochastic evolution equations in portfolio credit modelling. SIAM J. Finan. Math., 2, 627–664.

cappe
Cappé, O., Ryden, T, & Moulines, É. (2005). Inference
in Hidden Markov Models. Springer: New York.

chass_2019
Chassagneux, J., Szpruch, L., & Tse, A. (2019). Weak quantitative propagation of chaos via differential calculus on the space of measures. Ann. Appl. Probab., 32(3), 1929-1969.

delmoral
Del Moral, P. (2004).
Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications. Springer: New York.

delmoral1
Del Moral, P., Jacod, J. ,  & Protter, P. (2001). 
The Monte-Carlo method for filtering
with discrete-time observations.
Probab. Theory & Rel. Fields,  120, 346-368.

moral_2000
Del Moral, P., & Miclo, L. (2000). Branching and interacting particle systems. Approximations of Feynman-Kac formulae with applications to non-linear filtering. In Seminaire de Probabilites XXXIV, 1729, 1-145. 

dob
Dobramysl, U.,  Rudiger, S., Erban, R. (2016).  Particle-based multiscale modeling of calcium puff dynamics.  Multi. Mod. Sim.,  14, 997-1016.

dosreis
dos Reis, G., Smith, G. & Tankov, P. (2023).
Importance sampling for McKean-Vlasov SDEs.  Appl. Math. Comp., 453, 128078.

fearn
Fearnhead, P., Papaspiliopoulos,  O., & Roberts, G. O. (2008).  Particle filters for partially observed
diffusions. J. R. Stat. Soc.  Ser. B, 70, 755-777.

giles
Giles, M. B. (2008). Multilevel Monte Carlo path simulation. Op. Res., 56, 607-617.

giles1
Giles, M. B. (2015)  Multilevel Monte Carlo methods. Acta Numerica 24, 259-328.

haji
Haji-Ali,, A. L.,  Nobile, F.  &,  Tempone, R. (2016). Multi-index Monte Carlo: when
sparsity meets sampling.  Numer. Math.,  132, 767-806.

hein
Heinrich, S. (2001).
Multilevel Monte Carlo methods.
In Large-Scale Scientific Computing, (eds. S. Margenov, J. Wasniewski &
P. Yalamov), Springer: Berlin.







mlpf
Jasra, A., Kamatani, K., Law K. J. H. & Zhou, Y. (2017). 
Multilevel particle filters. SIAM J. Numer. Anal., 55, 3068-3096.

mlpf1
Jasra, A., Kamatani, K., Osei, P. P., & Zhou, Y. (2018). Multilevel particle filters: Normalizing Constant Estimation. Stat. Comp.,
28, 47-60.

ml_rev
Jasra, A.,  Law K. J. H.  & Suciu, C. (2020).  Advanced multilevel Monte Carlo.
Intl. Stat. Rev., 88,  548-579.

ubpf
Jasra, A.,  Law K. J. H.  & Yu, F. (2022).  Unbiased filtering of a class of partially observed diffusions.
Adv. Appl. Probab., 54, 661-687.

anti_mlpf
Jasra, A.,  Maama, M. & Ombao, H. (2024).  Antithetic multilevel particle filters.  
Adv. Appl. Probab. (to appear).


coup_clt
Jasra, A., & Yu, F. (2020). 
Central limit theorems for coupled particle filters.  Adv. Appl. Probab. 52, 942-1001.

mlpf_cts
Jasra,  A.,  Yu, F., & Heng, J. (2020).Multilevel particle filters for the non-linear filtering problem in continuous time. Stat. Comp., 30, 1381-1402.

koh
Kohatsu-Higa, A., & Ogawa, S. (1997). Weak rate of convergence for an Euler scheme of nonlinear SDE's.Monte Carlo Meth. Appl., 3(4), 327-345.


basic_method
Sznitman,  A. S. (1991). Topics in propagation of chaos. In Ecole d'Ete de Probabilites de Saint-Flour XIX.  1464, 165–251.

szp
Szpruch, L., Tan, S. & Tse, A. (2019). Iterative multilevel particle approximation of McKean-Vlasov SDEs.
Ann. Appl. Probab., 29, 2230-2265.


















