text:
\section*{Abstract}
In this paper, we study averaging principle for distribution dependent stochastic differential equations with localized $L^p$ drift. Using Zvonkin's transformation and estimates for Kolmogorov equations, we prove that the solution of the original system strongly and weakly converges to the solution of averaged system as the time scale $\varepsilon$ goes to zero. Moreover, we show the rates of strong and weak convergence depend on the singularity of the drift respectively.

\tableofcontents

\section{Introduction}
Consider the following distribution dependent stochastic differential equation (in short, DDSDE):
\begin{align}\label{in:NDDSDE}
{\mathord{{\rm d}}} X_t=b(t,X_t,\mu_t){\mathord{{\rm d}}} t+\Sigma(t,X_t,\mu_t){\mathord{{\rm d}}} W_t,
\end{align}
where $b: {\mathbb R}_+\times{\mathbb R}^d\times{\mathcal P}({\mathbb R}^d)\to{\mathbb R}^d$ and $\Sigma:{\mathbb R}_+\times{\mathbb R}^d\times{\mathcal P}({\mathbb R}^d)\to{\mathbb R}^d\otimes{\mathbb R}^d$ are measurable functions
and $\mu_t:={\mathcal L}(X_t)$ is the time marginal distribution or law of $X_t$. Here we write ${\mathcal P}({\mathbb R}^d)$ to mean the space of all probability measures over ${\mathbb R}^d$ and $(W_t)_{t\geqslant 0}$ is a $d$-dimensional standard Brownian motion on some stochastic basis $(\Omega,{\mathscr F},{\mathbf P}; ({\mathscr F}_t)_{t\geqslant 0})$. The infinitesimal generator of \eqref{in:NDDSDE} is given by
\begin{equation*}
{\mathscr L}_t\varphi(x)=a_{ij}(t,x,\mu_t)\partial_i\partial_j\varphi(x)+b(t,x,\mu_t)\cdot\nabla\varphi(x),
\end{equation*}
where $a_{ij}:=\frac{1}{2}\sigma_{ik}\sigma_{jk}$. Here and below we use the usual Einstein's
convention for summation: A repeated index in a product will be summed automatically.

Note that DDSDE \eqref{in:NDDSDE} is also called mean-field SDE or Mckean-Vlasov SDE in some literatures
if
\begin{equation}\label{MVb}
b(t,x,\mu):=\int_{\R^d}\hat{b}(t,x,y)\mu({\mathord{{\rm d}}} y),\quad \Sigma(t,x,\mu):=\int_{\R^d}\hat{\Sigma}(t,x,y)\mu({\mathord{{\rm d}}} y)
\end{equation}
for all $(t,x,\mu)\in\R_+\times\R^d\times{\mathcal P}(\R^d)$, where
$\hat{b}:{\mathbb R}_+\times{\mathbb R}^{2d}\to{\mathbb R}^d$ and $\hat{\sigma}:{\mathbb R}_+\times{\mathbb R}^{2d}\to{\mathbb R}^d\otimes{\mathbb R}^d$ are measurable functions.
It naturally appears in the studies on limiting behaviors of interacting particle systems and mean-field games. Roughly speaking,
DDSDE (or Mckean-Vlasov SDE) describes stochastic system whose evolution is determined not only by the microcosmic location of the particle, but also by the macrocosmic distribution. So it has attracted wide attention (see \cite{CD13,CD18b,CGPS20,Kac,Mc,Sz,Va}, in particular, \cite{CD21b} and references therein).

When $b$ and $\Sigma$ satisfy some continuous assumptions, there are numerous results devoted to studying the existence and uniqueness of this type of DDSDE \eqref{in:NDDSDE} (see \cite{Fu, HSS18, Wang2018} for examples). For the case that $b$ is only measurable, of at most linear growth and Lipschitz continuous with respect to $\mu$, when $\Sigma$ does not depend on $\mu$, is uniformly non-degenerate and Lipshitz continuous, by using the classical Krylov estimates, Mishura and Veretenikov \cite{MV2020} obtained the strong well-posedness of \eqref{in:NDDSDE}. When $\Sigma\equiv{\mathbb I}_{d\times d}$, the strong well-posedness were extended to local $L^q_tL^p_x$ drift by R\"ockner and Zhang in \cite{RZ21}. Moreover, by the entropy method and Girsanov's theorem, Lacker \cite{Lack2021} also obtained some well-posedness results for linear growth cases (see also \cite{La18}). Similar result was extended by Han in \cite{Ha22} for $L^q_tL^p_x$ drift. Notice that all the above results require that $b$ and $\Sigma$ depend linearly on the measure $\mu$. Namely, $b$ and $\Sigma$ satisfy \eqref{MVb}. In \cite{Zhao2020}, by some heat kernel estimates and Schauder-Tychonoff fixed point theorem, Zhao established the well-posedness results for DDSDE in general case \eqref{in:NDDSDE} when $b$ satisfies some $L^q_tL^p_x$ condition and the derivatives of $\Sigma$ with respect to $\mu$ are H\"older continuous. It should be noted that by Zvonkin's transformation of \cite{Zv} (see also \cite{Xi-Xi-Zh-Zh} and \cite{Zh-Zh1}) and entropy method, the last two authors together with Zhang in \cite{HRZ22} obtained the strong well-posedness for \eqref{in:NDDSDE} when $\Sigma$ does not dependent on $\mu$ and $b$ is in mixed $L^q_tL^p_x$ space.

In this work, we are interested in investigating averaging principle of the following
\iffalse
slow-fast system:
\begin{align*}
\begin{cases}
{\mathord{{\rm d}}} X_t^\varepsilon=b(Y^\varepsilon_t,X^\varepsilon_t,\mu^\varepsilon_t){\mathord{{\rm d}}} t+\sigma(X^\varepsilon_t){\mathord{{\rm d}}} W_t,\quad X^\varepsilon_0=\xi\in{\mathscr F}_0,\\
{\mathord{{\rm d}}} Y_t^\varepsilon=\frac{1}{\varepsilon}{\mathord{{\rm d}}} t,\quad Y^\varepsilon_0=0,
\end{cases}
\end{align*}
\fi
DDSDE with highly oscillating component
\begin{align}\label{MV1}
{\mathord{{\rm d}}} X_t^\varepsilon=b\left(\frac{t}{\varepsilon},X_t^\varepsilon,\mu_t^\varepsilon\right){\mathord{{\rm d}}} t+\sigma(X_t^\varepsilon){\mathord{{\rm d}}} W_t, \ \ X_0^\varepsilon=\xi\in{\mathscr F}_0,
\end{align}
where $\sigma:{\mathbb R}^d\to{\mathbb R}^d\otimes{\mathbb R}^d$ is a measurable function, $\mu^\varepsilon_t:={\mathcal L}(X^\varepsilon_t)$ is the distribution or law of $X^\varepsilon_t$, and small parameter $0<\varepsilon \ll1$ is the time scale.

Usually solving the original system \eqref{MV1} is
relatively difficult because of the highly oscillating component. Therefore, it is desirable to find a simplified system which simulates and predicts the evolution of the original system over a long time scale.
As well known, the highly oscillating component can be ``averaged'' out to produce such a simplified system under some suitable conditions, which is called
averaging principle.

More exactly, let
\begin{align}\label{in:KBM}
\bar{b}(x,\mu):=\lim\limits_{T\to\infty}\frac{1}{T}\int_0^T b(t,x,\mu){\mathord{{\rm d}}} t.
\end{align}
When $b$ is independent of $\mu$, $b$ is called
a {\em KBM-vector feld} (KBM stands for Krylov, Bogolyubov and Mitropolsky) if the convergence \eqref{in:KBM} exists uniformly with respect to $x$ in any
bounded subset on $\R^d$ (see e.g. \cite{SV1985}).
The averaging principle states that the solution of the original equation \eqref{MV1} converges,
as the time scale $\varepsilon$ goes to zero,
to that of the following averaged equation on any finite time intervals
\begin{align}\label{MV2}
{\mathord{{\rm d}}} X_t=\bar{b}(X_t,\mu_t){\mathord{{\rm d}}} t+\sigma(X_t){\mathord{{\rm d}}} W_t,\ \ X_0=\xi,
\end{align}
where $\mu_t$ stands for the distribution of $X_t$.

Averaging principle was firstly established for deterministic systems by Krylov, Bogolyubov and
Mitropolsky \cite{BM1961,KB1943}. Then it was extended to stochastic differential equations by
Khasminskii \cite{Khas1968}. After that,
extensive works concerning averaging principle for finite and infinite dimensional stochastic
differential equations were conducted; see e.g.
\cite{BK2004, Cerr2009, CF2009, CL2017, CL2021, DW2014, FW2012,
ML2020, Kif2004, MSV1991, PIX2021, RX2021, Skor1989, Vere1990, Vere1999,
WR2012} and the references therein.

\iffalse
Note that DDSDE \eqref{MV1} is also called mean-field SDE or Mckean-Vlasov SDE in some literatures
if
\begin{equation}\label{MVb}
b(t,x,\mu):=\int_{\R^d}b'(t,x,y)\mu({\mathord{{\rm d}}} y)
\end{equation}
for all $(t,x,\mu)\in\R_+\times\R^d\times{\mathcal P}(\R^d)$, where
$b':{\mathbb R}_+\times{\mathbb R}^{2d}\to{\mathbb R}^d$ is a measurable function.
It naturally appears in the studies on limiting behaviors of interacting particle systems and mean-field games. Roughly speaking,
DDSDE (or Mckean-Vlasov SDE) describes stochastic system whose evolution is determined not only by the microcosmic location of the particle, but also by the macrocosmic distribution. So it has attracted wide attention.

Under some monotone conditions, Wang \cite{Wang2018} proved the strong well-posedness and some functional inequalities of solutions for DDSDE.
When the coefficient is in the form of \eqref{MVb},
Mishura and Veretennikov \cite{MV2020} proved the existence and uniqueness of weak and strong solutions for Mckean-Vlasov SDE.
After that,
the results were extended; see e.g. \cite{Lack2021}. Moreover, R\"ockner and Zhang \cite{RZ21} obtained the well-posedness for
Mckean-Vlasov SDE with singular drifts in some localized $L_t^qL_x^p$ spaces. Recently, Zhao extended it to a more general case in \cite{Zhao2020}.
\fi

Recall that the strong convergence rate of averaging principle for slow-fast Mckean-Vlasov SDE
was established by the techniques of time discretization and Poisson equation in \cite{RSX2021}.
And as discussed in \cite{HLL2022}, the strong convergence rate of averaging principle for slow-fast Mckean-Vlasov SPDE was studied based on the variational approach and the technique of time discretization.
Note that the coefficients of the slow equation with fast variables were global Lipschitz continuous
with respect to slow variable in the above results.
Recently, the strong convergence without a rate for DDSDE with highly oscillating component
driven by fractional Brownian motion and standard Brownian motion was shown in \cite{SXW2022},
which requires that the drift term is continuous in slow variable.

Obviously, Lipschitz or continuous assumptions on $b$ are too strong for some applications.
There are a lot of interesting physical models only have bounded measurable or even singular $L^p$ interaction kernels $b$ and $b'$. For example, the rank-based interaction studied in \cite{La18} has an indicator kernel, which is discontinuous; the Biot-Savart law appearing in the vortex description of 2-dimensional incompressible Navier-Stokes equations has a singular kernel like $x^\perp/|x|^2$ (see e.g. \cite{Zh2020}). However, to the best of our knowledge, there is no result concerning the averaging principle both of DDSDE and SDE with $L^p$ drift.

From the above motivations, we consider the strong and weak convergence of averaging principle for DDSDE with $L^p$ drift in the present paper.
Moreover, we obtain the rate of strong and weak convergence, which is important for functional limit theorems in probability and homogenization in PDEs.
Throughout this paper we need the following conditions.

\begin{enumerate}[{\bf (H$_b^1$)}]
\item
Let $p_0\in(d\vee2,\infty)$ and assume that there is a nonnegative constant $\kappa_0$ such that for all $t\ge0$
and $\mu,\nu\in{\mathcal P}(\R^d)$
$$
|\mspace{-3mu}|\mspace{-3mu}| b(t,\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}+
\frac{|\mspace{-3mu}|\mspace{-3mu}| b(t,\cdot,\mu)-b(t,\cdot,\nu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}}{\|\mu-\nu\|_{var}}
\leqslant\kappa_0,$$
where
$\|\mu-\nu\|_{var}:=\sup\limits_{A\in{\mathscr B}({\mathbb R}^d)}|\mu(A)-\nu(A)|$ is total variation.\end{enumerate}

\begin{enumerate}[{\bf (H$_b^2$)}]
\item
There are functions
$\bar{b}:{\mathbb R}^{d}\times{\mathcal P}({\mathbb R}^d)\to{\mathbb R}^{d}$,
$\omega:{\mathbb R}_+\to{\mathbb R}_+$ and $H:{\mathbb R}^{d}\times{\mathcal P}({\mathbb R}^d)\to{\mathbb R}_+$ such that for all $(T,t,x,\mu)\in{\mathbb R}_+^2\times{\mathbb R}^{d}\times{\mathcal P}({\mathbb R}^d)$
\begin{align}\label{in:00}
\left|\frac{1}{T}\int_{t}^{T+t}(b(s,x,\mu)-\bar{b}(x,\mu)){\mathord{{\rm d}}} s\right|\leqslant \omega(T)H(x,\mu),
\end{align}
where $\lim\limits_{t\to\infty}\omega(t)=0$ and $\sup\limits_\mu|\mspace{-3mu}|\mspace{-3mu}| H(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}<\kappa_0$.
Here $p_0$ and $\kappa_0$ are as in {\bf(H$_b^1$)}.

\end{enumerate}
\begin{enumerate}[{\bf (H$_\sigma$)}]
\item
There are constants $p>d$, $\kappa_1>1$ and $\beta\in(0,1)$ such that for all $x,y,\xi\in{\mathbb R}^d$,
\begin{align*}
\kappa_1^{-1}|\xi|\leqslant |\sigma(x)\xi|\leqslant \kappa_1|\xi|,\quad |\mspace{-3mu}|\mspace{-3mu}| \nabla\sigma|\mspace{-3mu}|\mspace{-3mu}|_{p}\leqslant \kappa_1,
\end{align*}
and
$$
\|\sigma(x)-\sigma(y)\|_{HS}\leqslant \kappa_1|x-y|^{\beta},
$$
where $\|\cdot\|_{HS}$ is the Hilbert-Schmidt norm.
\end{enumerate}
Please see the definition of localized $L^p$ norm $|\mspace{-3mu}|\mspace{-3mu}|\cdot|\mspace{-3mu}|\mspace{-3mu}|_p$ in Section \ref{locBp}.
\begin{remark}
Notice that
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| \bar{b}(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}&\leqslant \frac{1}{T}\int_0^T|\mspace{-3mu}|\mspace{-3mu}| b(s,\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} s+|\mspace{-3mu}|\mspace{-3mu}| \frac{1}{T}\int_0^T\(b(s,\cdot,\mu)-\bar{b}(\cdot,\mu){\mathord{{\rm d}}} s\)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\\
&\leqslant \kappa_0+\omega(T)|\mspace{-3mu}|\mspace{-3mu}| H(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0},
\end{align*}
provided by conditions {\bf (H$^1_b$)} and {\bf (H$^2_b$)}. Taking $T\to\infty$, we have
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| \bar{b}(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\leqslant \kappa_0.
\end{align*}
Similarly, we have
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| \bar{b}(\cdot,\mu)-\bar{b}(\cdot,\nu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\leqslant \kappa_0\|\mu-\nu\|_{var}.
\end{align*}
Thus, the coefficient $\bar{b}$ satisfies the condition {\bf (H$^1_b$)} with the some constant $\kappa_0$.
\end{remark}
Under assumptions {\bf (H$_b^1$)} and {\bf (H$_\sigma$)}, for any initial value $\xi\in{\mathscr F}_0$,
it is well-known that there is a unique strong solution to DDSDE \eqref{MV1} (respectively, \eqref{MV2} ); see \cite{Zhao2020}\ and \cite{HRZ22}. The aim of this work is to show the following strong and weak convergence of averaging principle for DDSDE and SDE with $L^p$ drift.
\begin{theorem}\label{in:Main}
Under {\bf (H$_b^1$)}, {\bf (H$_b^2$)} and  {\bf(H$_\sigma$)},
 for any $T>0$ and $\ell\in(0,1)$,
there is a constant $C$, depending only on $\kappa_0,\kappa_1,T,d,\beta,p_0,p,\ell$, such that for any $\varepsilon>0$
\begin{align}\label{in:NMain}
\sup_{t\in[0,T]}\|\mu^\varepsilon_t-\mu_t\|_{var}\leqslant
 C \inf_{h>0}\Big(h^{\frac12-\frac{d}{2p_0}}
+\omega\left(\frac{h}{\varepsilon}\right)\Big)
\end{align}
and
\begin{align*}
{\mathbb E}\left(\sup_{t\in[0,T]}|X^\varepsilon_t-X_t|^{2\ell}\right)
\leqslant
 C\inf_{h>0}\left(\left(\omega(h/\varepsilon)\right)^2+h^{1-\frac{d}{p_0}}\right)^\ell.
 \end{align*}

\end{theorem}
\iffalse
\begin{remark}\rm
By Sobolev embedding \eqref{Eml00}, under the condition {\bf(H$_\sigma$)}, we have
\begin{align*}
\sigma\in {\mathbf C}^{1-\frac{d}{p_0}}.
\end{align*}
Hence, we can assume a priori that $\beta\geqslant 1-\frac{d}{p_0}$.

\end{remark}
\fi
When drift $b$ is independent of distribution, we have the following results, where the convergence rate is independent of $p_0$ and $\beta$.
\begin{theorem}\label{in:Main2}
Assume that
$$
b(t,x,\mu)\equiv b(t,x).
$$
 Under {\bf (H$_b^1$)}, {\bf (H$_b^2$)} and {\bf(H$_\sigma$)}, for any $T>0$, $\delta>0$ and $\ell\in(0,1)$,
there is a constant $C$, depending only on $\kappa_0,\kappa_1,T,d,p_0,p,\delta,\ell$, such that for any $\varepsilon>0$
$$
{\mathbb E}\left(\sup_{t\in[0,T]}|X^\varepsilon_t-X_t|^{2\ell}\right)
\leqslant C\inf_{h>0}\left(\left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2+h^{1-\delta}\right)^\ell.$$
\end{theorem}
\begin{remark}\rm
\begin{itemize}
\item[(i)] Since we use the Zvonkin transformation by parabolic equation,
when $\bar{b}=\bar{b}(t)=\bar{b}(\cdot,\mu_t)$ depends on time variable $t$,
the time regularity for solutions to this parabolic equation affects
the convergence rate (see \eqref{GG051} and Lemma \ref{LemEE01} for more details).
When $\bar{b}$ is independent of time, we can construct the Zvonkin transformation by elliptic equation. Hence, the convergence rates in Theorem \ref{in:Main} and \ref{in:Main2} are different. \\

\item[(ii)] Noting that $|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\lesssim\| f\|_\infty$ for all $p_0\in(1,\infty)$, all results in our paper are valid for $p_0=\infty$, in which the rate of convergence is
$$
\inf_{h>0}\left(\left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2+h^{1-\delta}\right)
$$
for any $\delta>0$. In particular, we obtain the convergence rate $\varepsilon^{\frac23-\delta}$ for a large number of examples (see e.g. Example \ref{Ex2} below), which is the same rate in \cite{RSX2021} when $\delta=0$. Moreover, the rate in \cite{HLL2022} is $\varepsilon^{\frac16}$.

\end{itemize}
\end{remark}

\begin{remark}\rm
These averaging principle results are also applicable to the following system
\begin{equation}\label{inseq}
{\mathord{{\rm d}}} X_t=\varepsilon b\left(t,X_t,{\mathcal L}(X_t)\right){\mathord{{\rm d}}} t+\sqrt{\varepsilon}\sigma(X_t){\mathord{{\rm d}}} W_t.
\end{equation}
With the help of time scaling $t\mapsto \frac{t}{\varepsilon}$, define $Y^\varepsilon_t:=X_{t/\varepsilon}$ and
$W^\varepsilon_t:=\sqrt{\varepsilon}W_{t/\varepsilon}$ for all $t\in\R_+$, we rewrite \eqref{inseq} as
\begin{equation*}
{\mathord{{\rm d}}} Y^\varepsilon_t=b\left(t/\varepsilon,Y^\varepsilon_t,{\mathcal L}(Y^\varepsilon_t)\right){\mathord{{\rm d}}} t+\sigma(Y^\varepsilon_t){\mathord{{\rm d}}} W^\varepsilon_{t}.
\end{equation*}
Then we can consider the following system
\begin{equation*}
{\mathord{{\rm d}}} \tilde{X}^\varepsilon_t=b\left(t/\varepsilon,\tilde{X}^\varepsilon_t,{\mathcal L}(\tilde{X}^\varepsilon_t)\right){\mathord{{\rm d}}} t+\sigma(\tilde{X}^\varepsilon_t){\mathord{{\rm d}}} W_{t}.
\end{equation*}
\end{remark}

Note that the drift of DDSDE and SDE in this paper is localized $L^{p_0}_x$ integrable,
we cannot use Gronwall's lemma or generalized Gronwall's lemma directly to
prove the convergence of $X^\varepsilon$ to $X$ as in \cite{HLL2022, SXW2022}. On the other hand, since the fast process $\{\frac{t}{\varepsilon}\}_{t\geq0}$ does not have the second elliptic generator, we cannot use the technique of the Poisson equation as in \cite{RSX2021}.
To overcome these difficulties, we use Zvonkin's transformation to remove the drift $b$
and employ the classical technique of time discretization.

More precisely, consider the following backward PDE related to \eqref{MV2}
\begin{align*}
\partial_tu+a_{ij}\partial_i\partial_ju-\lambda u+B\cdot\nabla u=B,\quad u(T)=0,
\end{align*}
where $B(t,x):=\bar{b}(x,\mu_t)$, $\lambda\ge0$ is the dissipative term.
Under {\bf (H$_\sigma$)} and {\bf (H$^1_b$)}, Xia et al \cite{Xi-Xi-Zh-Zh}
proved that for a sufficiently large number $\lambda$, there is a solution $u$ such that
\begin{align*}
|\nabla u(t,x)|\leqslant\frac{1}{2},\quad t\in[0,T], x\in{\mathbb R}^d.
\end{align*}
Hence, if we define $\Phi_t(x):=x+u(t,x)$, then $x\to\Phi_t(x)$ is a $C^1$ diffeomorphism of ${\mathbb R}^d$. Using It\^o's formula, $Y^\varepsilon_t:=\Phi_t(X^\varepsilon_t)$ and $Y_t:=\Phi_t(X_t)$ solve the following new SDEs:
\begin{align*}
{\mathord{{\rm d}}} Y^{\varepsilon}_t=
&
\lambda u(t,\Phi^{-1}_t(Y^\varepsilon_t)){\mathord{{\rm d}}} t+(\sigma^*\nabla\Phi_t(\Phi^{-1}_t(Y^\varepsilon_t))){\mathord{{\rm d}}} W_t\\
&
+\(b(t/\varepsilon,X^\varepsilon_t,\mu^\varepsilon_t)-\bar {b}(X^\varepsilon_t,\mu_t)\)\cdot\nabla \Phi_t\(X^{\varepsilon}_t\){\mathord{{\rm d}}} t
\end{align*}
and
\begin{align*}
{\mathord{{\rm d}}} Y_t=\lambda u(t,\Phi^{-1}_t(Y_t)){\mathord{{\rm d}}} t+(\sigma^*\nabla\Phi_t)(\Phi^{-1}_t(Y_t)){\mathord{{\rm d}}} W_t,
\end{align*}
where $\sigma^*$ is the transpose of $\sigma$ and $\Phi^{-1}_t$ is the inverse of $x\to\Phi_t(x)$.
Since these new systems have differentiable diffusion coefficients and the drifts are Lipshitz continuous,
we can use stochastic Gronwall's inequality. Please see the complete and rigorous formulation in Section \ref{Section5}.

The rest is how to use technique of time discretization to estimate the following crucial term
\begin{align}\label{in:NCD00}
{\mathbb E}\Big[\sup_{t\in[0,T]}\Big|\int_0^t\nabla \Phi_s\(X^{\varepsilon}_s\)\(b(\frac{s}{\varepsilon},X^\varepsilon_s,\mu^\varepsilon_s)-\bar{b}(X^\varepsilon_s,\mu_s)\){\mathord{{\rm d}}} s\Big|^2\Big].
\end{align}
In particular, we need to estimate
\begin{align}\label{in:NCC01}
\|\mu_t-\mu^\varepsilon_t\|_{var}
\end{align}
and the following difference for $B(t)=\nabla \Phi_tb(s/\varepsilon)$ and $B(t)=\nabla \Phi_t\bar{b}$:
\begin{align}\label{in:NCC00}
{\mathbb E}\Big|\int_0^t\(B(s,X^\varepsilon_{s},\mu^\varepsilon_{s})-B(s,X^\varepsilon_{\pi_h(s)},\mu^\varepsilon_{\pi_h(s)})\){\mathord{{\rm d}}} s\Big|^2,
\end{align}
where $\pi_h(s)$ is defined by $\pi_h(s)=s$ for $s\in[0,h)$ and
\begin{align*}
\pi_h(s):=kh,\quad s\in[kh,(k+1)h),~~\forall k\in{\mathbb N}.
\end{align*}
Since $b$ is not continuous, it is a challenging problem to give an estimate for \eqref{in:NCC00}.
Thanks to the important observation in \cite{DG20}, we deal with the problem by considering the time regularity for semi-group
generated by solutions to \eqref{MV1}. Basically, if $X^\varepsilon_t=W_t$ is a standard Brownian motion and $f$ is a discontinuous function, we have
\begin{align*}
&\quad{\mathbb E}\Big|\int_0^t\(f(W_{s})-f(W_{\pi_h(s)})\){\mathord{{\rm d}}} s\Big|^2\\
&=2{\mathbb E}\int_0^t\int_s^t\(f(W_{s})-f(W_{\pi_h(s)})\)\(f(W_{r})-f(W_{\pi_h(r)})\){\mathord{{\rm d}}} r{\mathord{{\rm d}}} s\\
&=2{\mathbb E}\int_0^t\(f(W_{s})-f(W_{\pi_h(s)})\)\Big({\mathbb E}^{{\mathscr F}_s}\int_s^t\(f(W_{r})-f(W_{\pi_h(r)})\){\mathord{{\rm d}}} r\Big){\mathord{{\rm d}}} s,
\end{align*}
where ${\mathscr F}_s:=\sigma(W_r,r\leqslant s)$ and ${\mathbb E}^{{\mathscr F}_s}$ is the conditional expectation with respect to ${\mathscr F}_s$.
Notice that by Markov property,
\begin{align*}
&\quad{\mathbb E}^{{\mathscr F}_s}\int_s^t\(f(W_{r})-f(W_{\pi_h(r)})\){\mathord{{\rm d}}} r\\
&={\mathbb E}^{{\mathscr F}_s}\int_s^{s+h}\(f(W_{r})-f(W_{\pi_h(r)})\){\mathord{{\rm d}}} r+\int_{s+h}^t\(P_{r-s}f(W_{s})-P_{\pi_h(r)-s}f(W_{s})\){\mathord{{\rm d}}} r,
\end{align*}
where $P_tf(x):={\mathbb E} f(x+W_t)$ is the semi-group generated by Brownian motion.
Note that the first term above is an infinitesimal one with respect to $h\to0$ and $P_tf$ is smooth.
That is, we transfer the function from a non-continuous one $f$ to a smooth one $P_tf$ in the second term.
Please see Section \ref{Section3} below for more details.
We mention that similar estimates are obtained in \cite{LL21} by stochastic sewing techniques.

For the difference \eqref{in:NCC01}, we employ a method based on Kolmogorov equation which is also used in \cite{RSX19}. Then, by time discretization again, we calculate the difference \eqref{in:NCC01} and obtain \eqref{in:NMain} (see Section \ref{Section4}).

\iffalse
More specifically,
for Theorem \ref{in:Main2}, the coefficient is independent of distribution,
the critical point in the analysis for convergence of $X^\varepsilon$ to $X$ is Lemma \ref{LemEE01}, which is based on Krylov's estimates and Girsanov's theorem. With the help of Lemma \ref{LemEE01} and Zvonkin's transformation, we obtain averaging principle for SDE with singular drift. When we consider DDSDE (Theorem \ref{in:Main}), another main obstacle is establishing the weak convergence. In particular, we show that $\mu^\varepsilon$ converges, as $\varepsilon$ goes to zero, to $\mu$ based on the Cauchy problem to Kolmogorov equation generated by $X_t,t\geq0$; see Theorem \ref{S4:Main}.
\fi

Now we illustrate our theoretical results by some examples. Firstly, the following example shows that the function $\omega(t)$ can be in the form $t^{-\alpha}$ with any power $\alpha\in(0,1]$.
\iffalse
\begin{Examples}\label{Ex0}\rm
Let
\begin{align}\label{Ex0eq}
b(t,x,\mu):=\left(\sum_{k=1}^\infty\left( {\mathbf{1}}_{[k,k+1]}(t)a_k\sin (2k\pi t)\right)\right)g(x,\mu),
\end{align}
where $g$ satisfies $\bf (H_1)$ and $|a_k|\leqslant f(k-1)$. Here nonnegative decreasing function $f$ satisfies $\lim\limits_{T\rightarrow\infty}\frac{f(T)}{T}=0$. Then it can be verified that
$$
\bar{b}(x,\mu)=g(x,\mu), \quad
\omega(t)=\frac{f(t)}{t},\quad \text{and}\quad H(x,\mu)=g(x,\mu).
$$
Note that $f$ can be $t^{-\alpha},\alpha>0$ and $e^{-t}$.
\end{Examples}
\fi
It also can be used in some singular interaction systems.
\begin{Examples}\label{Ex1}\rm
Consider the following DDSDE in ${\mathbb R}^d$
\begin{align*
{\mathord{{\rm d}}} X_t^\varepsilon
&
=\left(
\Big[(1+t/\varepsilon)^{-\alpha_1}+1\Big]\int_{\R^d} \frac{X_t^\varepsilon-y} {|X_{t}^{\varepsilon}-y|^{\alpha_2}}\mu_t^\varepsilon({\mathord{{\rm d}}} y)\right){\mathord{{\rm d}}} t+{\mathord{{\rm d}}} W(t)\\\nonumber
&
=:b(t/\varepsilon,X_{t}^{\varepsilon},\mu_{t}^{\varepsilon}){\mathord{{\rm d}}} t+{\mathord{{\rm d}}} W_t,
\end{align*}
where $\alpha_1>0$, $1<\alpha_2<2\wedge(1+\frac{d}{2})$ and $\mu_t^\varepsilon$ is the distribution of $X_t^\varepsilon$. It is clear that
the averaged equation is
\begin{align*}
{\mathord{{\rm d}}} X_t
&
=\left(\int_{\R^d}\frac{X_t-y} {|X_t-y|^{\alpha_2}}\mu_t({\mathord{{\rm d}}} y)\right){\mathord{{\rm d}}} t+{\mathord{{\rm d}}} W(t)\\
&
=:\bar{b}(X_t,\mu_t){\mathord{{\rm d}}} t+{\mathord{{\rm d}}} W(t)
\end{align*}
and
$$
\left|\frac{1}{T}\int_{t}^{t+T}\left(b(s,x,\mu)
-\bar{b}(x,\mu)\right){\mathord{{\rm d}}} s\right|
\leqslant\omega\left(T\right)(1-\alpha_1)^{-1}
\int_{\R^d} \frac{|x-y|} {|x-y|^{\alpha_2}} \mu({\mathord{{\rm d}}} y)
$$
for all $(T,t,x,\mu)\in\R^2\times\R\times{\mathcal P}(\R^d)$, where
\begin{eqnarray*}
\omega(t)=
\begin{cases}
t^{-(\alpha_1\wedge1)} & \text{for}~ \alpha_1\in(0,1)\cup(1,\infty)\\
t^{-1}\ln t           & \text{for}~ \alpha_1=1.
\end{cases}
\end{eqnarray*}
Then we have for any $\delta>0$
$$
\sup_{t\in[0,T]}\|\mu_t^\varepsilon-\mu_t\|_{var}\leqslant
C\varepsilon^{\frac{\alpha(2-\alpha_2)}{2+2\alpha-\alpha_2}-\delta}
$$
and
$$
\left[{\mathbb E}\left(\sup_{t\in[0,T]}|X^\varepsilon_t-X_t|^{2\ell}\right)
\right]^{\frac{1}{\ell}}
\leqslant C\varepsilon^{\frac{4\alpha-2\alpha\alpha_2}
{2+2\alpha-\alpha_2}-\delta}
$$
for any $0<\ell<1$, where $\alpha=\alpha_1\wedge1$ for $\alpha_1\in(0,\infty)$.
\end{Examples}

Then we give a more general example where the function $\omega(t)\asymp t^{-1}$.
We write $\omega(t)\asymp t$ to mean that there exists a constsnt
$C$ such that $C^{-1}t^{-1}\leqslant\omega(t)\leqslant Ct^{-1}$,
and put the check of assumptions {\bf (H$_b^1$)}-{\bf (H$_b^2$)} in Appendix.
\begin{Examples}\label{Ex2}\rm
Let $p_0>d\vee2$.
Consider the following DDSDE
\begin{align}\label{Ex2eq}
{\mathord{{\rm d}}} X_t^\varepsilon
=\left[\int_{\R} F\left(\sin(\xi t/\varepsilon),\int_{\R^d}\phi(X_t^\varepsilon,y)\mu_t^\varepsilon({\mathord{{\rm d}}} y)\right)\nu({\mathord{{\rm d}}} \xi)\right]
{\mathord{{\rm d}}} t+{\mathord{{\rm d}}} W_t,
\end{align}
where $F:[-1,1]\times{\mathbb R}^m\to{\mathbb R}^d$ is measurable and satisfies for some constant $L_F>0$
\begin{align}\label{in:EX2}
|F(t,0)|\leqslant L_F, ~~|F(t,x)-F(t,y)|\leqslant L_F|x-y|\quad \text{for all $(t,x,y)\in [-1,1]\times{\mathbb R}^{2m}$},
\end{align}
${\nu}$ is some finite measure on $\R$ satisfying
$$
\int_{{\mathbb R}\backslash\{0\}}\frac{1}{\xi}\nu({\mathord{{\rm d}}}\xi)<\infty
$$
and $\phi :{\mathbb R}^d\to{\mathbb R}^m$ is measurable and satisfies
$$
\sup\limits_y|\mspace{-3mu}|\mspace{-3mu}|\phi(\cdot,y)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}<\infty.
$$
Based on Theorem \ref{in:Main} and lemma \ref{lem61}, we have
$$
\sup_{t\in[0,T]}\|\mu_t^\varepsilon-\mu_t\|_{var}\leqslant
C\varepsilon^{\frac{1}{3}-\frac{2d}{9p_0-3d}}
$$
and
\begin{equation}\label{Ex2con}
\left[{\mathbb E}\left(\sup_{t\in[0,T]}|X^\varepsilon_t-X_t|^{2\ell}\right)
\right]^{\frac{1}{\ell}}
\leqslant C\varepsilon^{\frac{2}{3}\left(1-\frac{2d}{3p_0-d}\right)}
\end{equation}
for any $0<\ell<1$ and $T>0$, where $X_t$ is the solution of
\begin{align*}
{\mathord{{\rm d}}} X_t
&
=\Big[\frac{1}{2\pi}\int_0^{2\pi} F\left(\sin\tau,\int_{\R^d}\phi(X_t,y)\mu_t({\mathord{{\rm d}}} y)\right)\nu(\R\setminus\{0\}){\mathord{{\rm d}}}\tau\\
&\qquad
+F\left(0,\int_{\R^d}\phi(X_t,y)\mu({\mathord{{\rm d}}} y)\right)\nu(\{0\})\Big]{\mathord{{\rm d}}} t+{\mathord{{\rm d}}} W(t).
\end{align*}
Here $\mu_t$ is the distribution of $X_t$ for all $t\geq0$.
\end{Examples}

\begin{remark}\rm
In this paper, we only consider the case when $\sigma$ does not depend on time $t$ and measure $\mu$. We hope that in future work
we can use our framework here to study the time-inhomogeneous case $\sigma=\sigma(t,x)$ or an even more general case $\sigma=\sigma(t,x,\mu)$ under the following standard condition of $\sigma$
\begin{enumerate}[{\bf (H$_\sigma^{'}$)}]
\item
There are functions $\omega_\sigma$ and $H_\sigma$ such that
\begin{align*
\frac{1}{T}\int_{t}^{T+t}\|\sigma(t,x,\mu)-\bar{\sigma}(x,\mu)\|_{HS}^2{\mathord{{\rm d}}} t\leqslant \omega_\sigma(T)H_\sigma(x,\mu)
\end{align*}
for all $(T,t,x)\in(\R_+^2\times\R^d)$,
where $\lim\limits_{t\to\infty}\omega_\sigma(t)=0$ and $\sup\limits_\mu|\mspace{-3mu}|\mspace{-3mu}| H_\sigma(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{\frac{p_0}{2}}<\infty$.
\end{enumerate}

\iffalse
then the semigroup $P_{s,t}^\sigma$ generated by $$Z_{s,t}^x=x+\int_s^t\sigma(r,Z_{s,r}^x){\mathord{{\rm d}}} W(r),~0\leqslant s\leqslant t$$ is time-inhomogeneous
(see Section \ref{sec:3.1} for details concerning homogeneous case).
Therefore, we cannot obtain the regularity of time for $P_{s,t}^\sigma$
by the method in Lemma \ref{lemmm03}.
To overcome this difficulty, we use the technique of Levi's freezing coefficients
and Duhamel's formula (see Appendix for details) to have
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| P^\sigma_{s,t_1}\varphi-P^\sigma_{s,t_2}\varphi|\mspace{-3mu}|\mspace{-3mu}|_q\lesssim [(t-s)^{\frac{\beta}{2}}s^{-\frac{\beta}{2}}]s^{\frac{d}{2q}-\frac{d}{2p}}.
\end{align*}
Now, consider the following DDSDE driven by multiplicative noise with highly oscillating component
\begin{align}\label{ZZ0311}
{\mathord{{\rm d}}} X^\varepsilon_t=b\left(\frac{t}{\varepsilon},X^\varepsilon_t,\mu^\varepsilon_t\right){\mathord{{\rm d}}} t+\sigma\left(\frac{t}{\varepsilon},X^\varepsilon_t\right){\mathord{{\rm d}}} W_t.
\end{align}
Assume that $\bar{\sigma}(x):=\lim\limits_{T\to\infty}\frac{1}{T}\int_0^T \sigma(t,x){\mathord{{\rm d}}} t$ for all $x\in\R^d$.
If, {\bf(H$_\sigma$)}, {\bf (H$_b^1$)}, {\bf (H$_b^2$)} and also {\bf (H$_\sigma^{'}$)} hold,
we can still use our method in this paper to obtain the convergence rate for \eqref{ZZ0311}.
Here
\begin{enumerate}[{\bf (H$_\sigma^{'}$)}]
\item
There are functions $\omega_\sigma$ and $H_\sigma$ such that
\begin{align*
\frac{1}{T}\int_{t}^{T+t}\|\sigma(t,x)-\bar{\sigma}(x)\|_{HS}^2{\mathord{{\rm d}}} t\leqslant \omega_\sigma(\frac{1}{T})H_\sigma(x)
\end{align*}
for all $(T,t,x)\in(\R_+^2\times\R^d)$,
where $\lim\limits_{t\to0}\omega_\sigma(t)=0$ and $\sup\limits_\mu|\mspace{-3mu}|\mspace{-3mu}| H_\sigma|\mspace{-3mu}|\mspace{-3mu}|_{\frac{p_0}{2}}<\infty$,
\end{enumerate}
In this case, the order of convergence rates for SDE and DDSDE are the same.
For simplicity, we only consider the case when $\sigma$ is time-independent in what follows.
\fi

\end{remark}

{\bf Structure of the paper.}

In Section \ref{Section2}, we first introduce the localized Bessel potential space
and some basic properties used in this paper. Then we show some well-posedness results
and priori estimates about parabolic and elliptic PDEs for later establishing
the time regularity of solutions in Section \ref{Section3} and
Zvonkin's transformation in Section \ref{Section5}.

In Section \ref{Section3}, we study the time regularity of solutions to
both SDE and PDE under the $L^p$ framework.
As usual, we use the technique of time discretization
to obtain the convergence for averaging principle.
That is we need to show a time discretization's type estimate
of \eqref{in:NCC00},
and then estimate \eqref{in:NCD00} by condition {\bf (H$^2_b$)} (see Section \ref{Sec3.1}).
For \eqref{in:NCC00}, we first asume that $Z_t^\varepsilon$ is a solution to
the SDE \eqref{MV1} without drift, that is to say $b\equiv0$,
by the semi-group property from condition {\bf (H$_\sigma$)}
and heat kernel estimates we estimate \eqref{in:NCC00}
for any localized $L^p$ integrable $B$ with $p>(d\vee2)$. Then it follows
from the Girsanov transform that \eqref{in:NCC00} holds for $X_t^\varepsilon$,
where $X_t^\varepsilon$ is the solution to \eqref{MV1} (see Section \ref{sec:3.1}).
Moreover, the difference $\|\mu_s-\mu_{\pi_h(s)}\|_{{\mathrm{var}}}$ comes from Duhamel's formula
and time regularity of solutions to SDEs in Section \ref{Sec3.1}.
We also establish the time regularity of the gradient of solutions for parabolic PDE in Section \ref{SubS32},
which is used in Section \ref{Section4} and Section \ref{Section5} respectively.

In Section \ref{Section4}, by a trick of observing the relation between the solution to SDE and
the solution to Kolmogorov equation,
we give the weak convergence rate of \eqref{in:NCC01},
where the key point is the estimates of time regularity obtained in Section \ref{SubS32}.

In Section \ref{Section5},
we give the proofs to our main results Theorem \ref{in:Main}
and \ref{in:Main2} by Zvonkin's
transformation.
\vspace{5mm}

{\bf Notations.}

Throughout this paper, we denote by $|\cdot|$ the Euclidean norm on ${\mathbb R}^N$
with dimension $N\in{\mathbb N}$ whose value may change in different places.
For $j\in{\mathbb N}\cup\{0\}$, we use the notation $\nabla^j$ to denote
the $j$th order derivative. Moreover, denote by $C^\infty_0$
the function spaces of all smooth functions with compact support.
We write $C^k_b$ (respectively, $C^\infty_b$) to mean the space of all smooth functions
with bounded $j$th derivatives for all integer $j\in[0,k]$ (respectively, $j\in{\mathbb N}\cup\{0\}$).
Let $[\alpha]$ be the largest integer which is smaller than $\alpha$
for any constant $\alpha>0$.
We use $C$ with or without subscripts to denote an unimportant constant,
whose value may change on different occasions.
We also denote by $:=$ a way of definition.
By $A\lesssim B$
we mean that for some unimportant constant $C\ge1$,
\begin{align*}
A\leqslant CB.
\end{align*}

\section{Preliminaries}\label{Section2}

In this section we show some auxiliary results for later use.
\subsection{Localized Bessel potential space}\label{locBp}
For any $(\alpha,p)\in\R\times(1,\infty)$, we write
$$
H^{\alpha,p}:=(\mathbb I-\Delta)^{-\alpha/2}
\left(L^p(\R^d)\right)
$$
to mean the usual Bessel potential space with the norm given by
$
\|f\|_{\alpha,p}:=\|(\mathbb I-\Delta)^{\alpha/2}f\|_p,
$
where $\|\cdot\|_p$ is the usual $L^p$-norm in $\R^d$.
Here $(\mathbb I-\Delta)^{\alpha/2}f$ is defined through
Fourier's transform
$$
(\mathbb I-\Delta)^{\alpha/2}f:={\mathcal F}^{-1}\left((1+|\cdot|^2)^{\alpha/2}
{\mathcal F} f\right).
$$
Note that for any $n\in\mathbb N} \def\kk{\kappa} \def\m{{\bf m}$ and $p\in(1,\infty)$, an equivalent norm in $H^{n,p}$
is given by (see e.g. \cite{St70})
$$
\|f\|_{n,p}=\|f\|_p+\|\nabla^nf\|_p.
$$

Let $\chi\in C_0^\infty(\R^d)$
such that $\chi(x)=1$ for $|x|\leq1$,  $\chi(x)=0$
for $|x|>2$ and $|\chi(x)|\leq1, \forall x\in\R^d$. Define
$$
\chi_r(x):=\chi(x/r), \quad \chi_r^z(x):=\chi_r(x-z)
$$
for all $r>0$ and $z\in\R^d$. Given $r>0$, we introduce the
following localized $H^{\alpha,p}$-space:
$$
\widetilde{H}^{\alpha,p}:=\left\{f\in H^{\alpha,p}_{\rm{loc}}(\R^d),|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\alpha,p}:=\sup_z\|\chi_r^zf\|_{\alpha,p}<\infty\right\}.
$$
When $\alpha=0$, we simply write
$$
|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p:=|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{0,p}.
$$

For $T>0$, $p,q\in(1,\infty)$ and $\alpha\in\R$, we set
\begin{align*}
  {\mathbb L}_q^p(T):=L^q\([0,T];L^p\),\quad
  {\mathbb H}_q^{\alpha,p}(T):=L^q\([0,T];H^{\alpha,p}\).
\end{align*}
Now we introduce the localized space $\widetilde{{\mathbb H}}_q^{\alpha,p}(T)$ with the following norm
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}_q^{\alpha,p}(T)}:=
\sup_{z\in\R^d}\|\chi_r^zf\|_{{\mathbb H}_q^{\alpha,p}(T)}<\infty.
\end{align*}
By a finite covering technique, it can be verified  that the definition of $\widetilde{{\mathbb H}}_q^{\alpha,p}$
does not depend on the choice of $r$ (see \cite[Section 2]{Xi-Xi-Zh-Zh}).
Note that
\begin{align*}
  L^q\([0,T];\widetilde{H}^{\alpha,p}\)\subset
  \widetilde{{\mathbb H}}_q^{\alpha,p}(T).
\end{align*}
For simplicity, we define
\begin{align*}
\widetilde{{\mathbb H}}^{\alpha,p}(T):=
L^\infty([0,T];\widetilde{H}^{\alpha,p}),
\quad
\widetilde{{\mathbb L}}^p(T):=\widetilde{{\mathbb L}}^p_\infty(T),
\quad\text{and}\quad {\mathbb L}^\infty_T:=L^\infty([0,T]\times{\mathbb R}^d).
\end{align*}
It follows from H\"older's inequality that for any $1\leqslant p_2\leqslant p_1\leqslant\infty$
\begin{align*
L^{p_1}\subset\widetilde{L}^{p_1}\subset\widetilde{L}^{p_2}.
\end{align*}

Denote by $C^\alpha$ the H\"older space which consists of all functions $g$ for which the norm
\begin{align*}
\|g\|_{C^\alpha}:=\sum_{|\beta|\leqslant[\alpha]}
\|\nabla^\beta g\|_{L^\infty}
+\sum_{|\beta|=[\alpha]}
\left[\nabla^\beta g\right]_{C^{\alpha-[\alpha]}}
\end{align*}
is finite, where
\begin{align*}
  [g]_{C^{\alpha-[\alpha]}}:=\sup_{x,y\in\R^d \atop x\neq y}
  \frac{|g(x)-g(y)|}{|x-y|^{\alpha-[\alpha]}}.
\end{align*}

\begin{lemma}[Embedding lemma]\label{Eml}
Let $1<p<\infty$.
Then we have
\begin{align*
\widetilde{H}^{\alpha,p}\subset C^{\alpha-d/p}
\end{align*}
and
\begin{align*}
  \widetilde{{\mathbb H}}^{\alpha,p}\subset L^\infty\left([0,T];C^{\alpha-d/p}\right)
\end{align*}
provided $\alpha>d/p$.
\end{lemma}
\begin{proof}
It follows from Sobolev's embedding theorem that $H^{\alpha,p}\subset C^{\alpha-d/p}$ if $\alpha>d/p$.
Note that
\begin{align*}
\|g\|_{C^{\alpha-d/p}}\leqslant\sup_{z}\|\chi_r^zg\|_{C^{\alpha-d/p}}
\end{align*}
for all $g\in C^{\alpha-d/p}$ and $r>0$. Therefore, we have
\begin{align*}
\|g\|_{C^{\alpha-d/p}}
\leqslant\sup_{z}\|\chi_r^zg\|_{C^{\alpha-d/p}}
\lesssim \sup_{z}\|\chi_r^zg\|_{H^{\alpha,p}}
=|\mspace{-3mu}|\mspace{-3mu}| g|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{H}^{\alpha,p}}.
\end{align*}
Moreover, for all $f\in\widetilde{{\mathbb H}}^{\alpha,p}$ one sees that
\begin{align*}
\sup_{t\in[0,T]}\|f(t)\|_{C^{\alpha-d/p}}
\leqslant\sup_{t\in[0,T]}\sup_{z}\|\chi_r^zf(t)\|_{C^{\alpha-d/p}}
\lesssim \sup_{t\in[0,T]}\sup_{z}\|\chi_r^zf(t)\|_{H^{\alpha,p}}
=|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}^{\alpha,p}}.
\end{align*}
\end{proof}

Now we introduce the local Hardy-Littlewood maximal function, which is defined by
\begin{align*}
{\mathcal M} f(x):=\sup_{r\in(0,1)}\frac{1}{|B_r|}\int_{B_r}f(x+y){\mathord{{\rm d}}} y.
\end{align*}

The following result is taken from Lemma 2.1 in \cite{Xi-Xi-Zh-Zh} (see also \cite[Appendix A]{CD08}).

\begin{lemma}
\begin{enumerate}[(i)]
\item There is a constant $C=C(d)>0$, such that for any $f\in L^\infty({\mathbb R}^d)$ with $\nabla f\in L^1_{\rm{loc}}({\mathbb R}^d)$
and for Lebesgue-almost all $x,y\in{\mathbb R}^d$,
\begin{align}\label{Pre1:Mnf}
|f(x)-f(y)|\leqslant C|x-y|\({\mathcal M}|\nabla f|(x)+{\mathcal M}|\nabla f|(y)+\|f\|_\infty\).
\end{align}
\item For any $p\in(1,\infty)$, there is a constant $C=C(d,p)>0$ such that for all $f\in \widetilde{L}^p$,
\begin{align}\label{Pre1:Mix}
|\mspace{-3mu}|\mspace{-3mu}| {\mathcal M} f|\mspace{-3mu}|\mspace{-3mu}|_{p}\leqslant C|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}.
\end{align}
\end{enumerate}
\end{lemma}

\subsection{Parabolic equation}
In order to study DDSDE, we consider the following second order parabolic PDE in ${\mathbb R}_+\times{\mathbb R}^d$:
\begin{align}\label{Pre4:PDE}
\partial_tu=a_{ij}\partial_i\partial_ju-\lambda u+b\cdot\nabla u+f,\ \ u(0)=\varphi,
\end{align}
where $\lambda\ge0$, $a=(a_{ij}):\R_+\times\R^d\rightarrow\R^d\otimes\R^d$
is a symmetric matrix-valued Borel measurable function satisfying (H$_a$), i.e.
\begin{enumerate}
  \item[(H$_a$)] there exist constants $c_0>0$ and $\theta\in(0,1)$ such that
  \begin{equation*}
    c_0^{-1}|\xi|\leqslant|a(x)\xi|\leqslant c_0|\xi|,\quad
    \|a(x)-a(y)\|_{HS}\leqslant c_0|x-y|^\theta
  \end{equation*}
  for all $\xi\in\R^d$ and $x,y\in\R^d$,
\end{enumerate}
and
$b:\R\times\R^d\rightarrow\R^d$ is vector-valued Borel measurable function.
Firstly, we introduce the definition of the solution to PDE \eqref{Pre4:PDE}.
\begin{definition}\label{Pre4:Def}\rm
Let $T>0$, $p,q\in(1,\infty)$, $\lambda\ge0$,
$b,f\in \widetilde{{\mathbb L}}^p_q(T)$ and $\varphi\in C_b^\infty$.
We call a function $u$ with
$\partial_tu\in\widetilde{{\mathbb L}}_q^p(T)$ and
$u\in \widetilde{{\mathbb H}}^{2,p}_q(T)$ a {\em solution} of PDE \eqref{Pre4:PDE} if for Lebesgue almost all $(t,x)\in{\mathbb R}_+\times{\mathbb R}^d$,
\begin{align*}
u(t,x)=\int_0^t\(a_{ij}\partial_i\partial_ju(s,x)-\lambda u(s,x)+b\cdot\nabla u(s,x)+f(s,x)\){\mathord{{\rm d}}} s+\varphi(x).
\end{align*}
\end{definition}

\begin{remark}\rm
For any $\chi\in C_0^\infty(\R^d)$ and $f\in\widetilde{\mathbb H}_q^{\alpha,p}(T)$,
by the definition of localized spaces $\widetilde{\mathbb H}_q^{\alpha,p}(T)$,
we have $\chi f\in{\mathbb H}_q^{\alpha,p}(T)$. Hence for any solution $u$ of PDE
\eqref{Pre4:PDE} in the sense of Definition \ref{Pre4:Def}, $\chi u$ is H\"older continuous
on $[0,T]\times\R^d$
when $d/p+2/q<2$ provided by \cite[Lemma 10.2]{KR05}. Moreover,
$\nabla (\chi u)$ is H\"older continuous on $[0,T]\times\R^d$
if $d/p+2/q<1$. In view of the arbitrariness
of cut-off function $\chi$, $u$ (respectively, $\nabla u$) is H\"older continuous
on $[0,T]\times\R^d$ if $d/p+2/q<2$ (respectively, $d/p+2/q<1$).
\end{remark}

The following property of the solution $u$ comes from \cite[Theorem 3.2]{Xi-Xi-Zh-Zh}.
\begin{lemma}\label{lem210}
Let $T>0$, $p,q\in(1,\infty)$ with $2/q+d/p<1$,
$\lambda\ge0$, $b\in \widetilde{{\mathbb L}}^{p}_{q}(T)$. Set
$$
\Theta:=(d,T,p,q,\|b\|_{\widetilde{{\mathbb L}}^{p}_{q}(T)},c_0,\theta).
$$ Then there is a constant $\lambda_0=\lambda_0(\Theta)$ such that for all $\lambda\geqslant\lambda_0$ and $f\in \widetilde{{\mathbb L}}^{p}_{q}(T)$ and $\varphi\in \widetilde{H}^{2,p}$,
there is a unique solution $u$ to PDE \eqref{Pre4:PDE} on $[0,T]$ in the sense of Definition \ref{Pre4:Def} satisfying  that for any $\alpha\in[0,2)$, $p'\in[p,\infty]$, $q'\in[q,\infty]$ with
$$
\beta:=2-\alpha+\frac{2}{q'}+\frac{d}{p'}-\(\frac{2}{q}+\frac{d}{p}\)>0,
$$
there is a constant  $C=C(\Theta,\alpha,p',q')>0$ such that for all $\lambda\geqslant\lambda_0$
\begin{align}\label{CCA02}
\lambda^{\frac{\beta}{2}}
|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}^{\alpha,p'}_{q'}(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\partial_tu|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\nabla^2 u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}
\leqslant C\(|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}+|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p}\).
\end{align}
\end{lemma}
\begin{remark}\label{Re001}\rm
By Lemma \ref{Eml}, we have $u\in L^\infty([0,T]; {\mathbf C}^\gamma)$ for any $\gamma\in(1,2-2/q-d/p)$.
\end{remark}
\begin{proof}
Notice that $u$ is a solution to PDE \eqref{Pre4:PDE} with
$u(0)=\varphi$ in the sense of Definition \ref{Pre4:Def}
if and only if $\bar{u}:=u-\varphi$ is a solution to PDE
\eqref{Pre4:PDE} with $\bar{u}(0)=0$ and
$f=f+a_{ij}\partial_i\partial_j\varphi-\lambda\varphi-b\cdot\nabla \varphi$.
Based on Lemma \ref{Eml}, we have
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| f+a_{ij}\partial_i\partial_j\varphi-\lambda\varphi-b\cdot\nabla \varphi|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}&\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}+ |\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p}+ |\mspace{-3mu}|\mspace{-3mu}| b|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}\|\nabla\varphi\|_\infty\\
&\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}+(1+ |\mspace{-3mu}|\mspace{-3mu}| b|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)})|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p},
\end{align*}
and complete the proof by \cite[Theorem 3.2]{Xi-Xi-Zh-Zh}.
\end{proof}

With the help of prior estimate \eqref{CCA02}, we obtain
the well-posedness of PDE \eqref{Pre4:PDE} for any $\lambda\geq0$.
\begin{proposition}\label{well-PDE}
Let $T>0$, $p,q\in(1,\infty)$ with $2/q+d/p<2$, $\lambda\geq0$, $b\in\widetilde{{\mathbb L}}_q^p(T)$.
Then for all $f\in \widetilde{{\mathbb L}}^{p}_{q}(T)$ and $\varphi\in \widetilde{H}^{2,p}$
there is a unique solution $u$ to PDE \eqref{Pre4:PDE} on $[0,T]$
in the sense of Definition \ref{Pre4:Def} such that
\begin{equation}\label{CCA022}
|\mspace{-3mu}|\mspace{-3mu}|\nabla u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_T^\infty}
+|\mspace{-3mu}|\mspace{-3mu}|\partial_tu|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}
+|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}_q^{2,p}(T)}\leqslant
C\left(|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{H}^{2,p}}\right),
\end{equation}
where $C=C(\Theta,\lambda)$.
\end{proposition}
\begin{proof}
By standard continuity method, it suffices to show the prior estimate \eqref{CCA022} for \eqref{Pre4:PDE}. To this end, we rewrite \eqref{Pre4:PDE} to
\begin{equation*
\partial_tu=a_{ij}\partial_i\partial_ju-(\lambda+\lambda_0)u+b\cdot\nabla u
+f+\lambda_0u,
\end{equation*}
where $\lambda_0$ is as in Lemma \ref{lem210}. In view of Lemma \ref{lem210}, we have
\begin{align}\label{RePDE1}
(\lambda+\lambda_0)^{\frac{\beta}{2}}|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}^{\alpha,p}_{q'}}
+|\mspace{-3mu}|\mspace{-3mu}|\partial_tu|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\nabla^2u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}\leqslant
C\left(|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}
+\lambda_0|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p}\right),
\end{align}
where $C=C(\Theta,\alpha,q')>0$, $\beta=2-\alpha+\frac{2}{q'}-\frac{2}{q}>0$.
Taking $q'=\infty$ in \eqref{RePDE1}, then we have
\begin{align*
(\lambda+\lambda_0)^{\frac{\beta}{2}}\sup_{t\in[0,T]}
|\mspace{-3mu}|\mspace{-3mu}| u(t)|\mspace{-3mu}|\mspace{-3mu}|_p\leqslant
C\left(|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}
+\lambda_0\left(\int_0^T|\mspace{-3mu}|\mspace{-3mu}| u(t)|\mspace{-3mu}|\mspace{-3mu}|_p^q{\mathord{{\rm d}}} t\right)^{\frac{1}{q}}
+|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p}\right).
\end{align*}
It follows from Gronwall's Lemma that
\begin{align}\label{RePDE3}
\sup_{t\in[0,T]}|\mspace{-3mu}|\mspace{-3mu}| u(t)|\mspace{-3mu}|\mspace{-3mu}|_p\leqslant C\left(
|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p}\right),
\end{align}
where $C$ depends on $\Theta,\alpha,\lambda$. Combining
\eqref{RePDE1} and \eqref{RePDE3}, we obtain for $1<\alpha<2-2/q$
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}^{\alpha,p}}
+|\mspace{-3mu}|\mspace{-3mu}|\partial_tu|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}
+|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}_q^{2,p}(T)}\leqslant
C\left(|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}_q^p(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_{2,p}\right).
\end{align*}
and complete the proof by Lemma \ref{Eml}.
\end{proof}

\subsection{Elliptic equation}
Now we consider the following second order elliptic PDE in ${\mathbb R}^d$:
\begin{align}\label{Pre3:PDE}
a_{ij}\partial_i\partial_ju-\lambda u+b\cdot\nabla u=f,
\end{align}
where $\lambda\ge0$, $a=(a_{ij}):\R^d\rightarrow\R^d\otimes\R^d$
is a symmetric matrix-valued Borel measurable function satisfying (H$_a$) and
$b:\R^d\rightarrow\R^d$ is vector-valued Borel measurable function.
Firstly, we introduce the definition of the solution to PDE \eqref{Pre3:PDE}.
\begin{definition}\label{Pre3:Def}\rm
Let $p\in(1,\infty)$, $\lambda,T\ge0$ and $b,f\in \widetilde{L}^p$. We call
$u\in\widetilde{H}^{2,p}$ a {\em solution} of PDE \eqref{Pre3:PDE} if for Lebesgue almost all $x\in{\mathbb R}^d$,
\begin{align*}
a_{ij}\partial_i\partial_ju(x)-\lambda u(x)+b(x)\cdot\nabla u(x)=f(x).
\end{align*}
\end{definition}
As a corollary of Lemma \ref{lem210}, we have the following results (see \cite[Theorem 7.6]{Xi-Zh}).

\begin{lemma}
Assume $b\in \widetilde{L}^{p}$ with some $p>d$.
There are constants $\lambda_0=\lambda_0(d,|\mspace{-3mu}|\mspace{-3mu}| b|\mspace{-3mu}|\mspace{-3mu}|_p,p,c_0,\theta)$
and $C=C(d,|\mspace{-3mu}|\mspace{-3mu}| b|\mspace{-3mu}|\mspace{-3mu}|_p,p,p',c_0,\theta)$ such that for any $\lambda\geqslant\lambda_0$
and $f\in\widetilde{L}^{p}$, there exists a unique solution
$u$ to the PDE \eqref{Pre3:PDE} in the sense of Definition \ref{Pre3:Def} such that
\begin{align}\label{CC02}
\lambda^{\frac{\beta}{2}}
|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{H}^{\alpha,p'}}
+|\mspace{-3mu}|\mspace{-3mu}|\nabla^2 u|\mspace{-3mu}|\mspace{-3mu}|_{p}\leqslant C|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p},
\end{align}
where $\alpha\in[0,2),p\in[p,\infty], p'\in[p,\infty]$ with
$\beta:=2-\alpha+\frac{d}{p'}-\frac{d}{p}>0$.
\end{lemma}
\begin{proof}
As usual, it suffices to show a priori estimate \eqref{CC02}.
Let $T>0$, $u$ be a solution to \eqref{Pre3:PDE}
and $\phi$ be a nonnegative and nonzero smooth function on $[0,\infty)$ with $\phi(0)=0$.
Define $\tilde{u}(t,x):=\phi(t)u(x)$. Then, one sees that
$\tilde{u}$ is a solution to the following parabolic
equation in the sense of Definition \ref{Pre4:Def}:
\begin{align*}
\partial_t\tilde{u}=a_{ij}\partial_i\partial_j\tilde{u}-\lambda\tilde{u}+b\cdot\nabla\tilde{u}-\phi f+\phi'u,\quad \tilde{u}(0)=0.
\end{align*}
Based on \eqref{CCA02}, we have for any $\alpha\in[0,2),p'\in[p,\infty]$ with
$\beta:=2-\alpha+\frac{d}{p'}-\frac{d}{p}>0$,
\begin{align*}
\lambda^{\frac{\beta}{2}}
|\mspace{-3mu}|\mspace{-3mu}| \tilde{u}|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb H}}^{\alpha,p'}(T)}
+|\mspace{-3mu}|\mspace{-3mu}|\nabla^2 \tilde{u}|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}\leqslant C|\mspace{-3mu}|\mspace{-3mu}|\phi' u-\phi f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)},
\end{align*}
which implies that
\begin{align*}
\lambda^{\frac{\beta}{2}}
|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{H}^{\alpha,p'}}
+|\mspace{-3mu}|\mspace{-3mu}|\nabla^2 u|\mspace{-3mu}|\mspace{-3mu}|_{p}\leqslant C\|\phi\|^{-1}_\infty
\left(\|\phi'\|_\infty |\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_p+\|\phi\|_\infty |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p\right),
\end{align*}
where $\|\phi\|_\infty:=\sup_{t\in[0,T]}|\phi(t)|$. Noting that $|\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_p\leqslant |\mspace{-3mu}|\mspace{-3mu}| u|\mspace{-3mu}|\mspace{-3mu}|_{p'}$,
we obtain \eqref{CC02} and complete the proof.

\end{proof}

\section{Analysis of time regularity}\label{Section3}
{ In this section,
 letting $T>0$, we assume that
 \begin{align}\label{NBB03}
 B\in \widetilde{{\mathbb L}}^{p_0}(T)\text{ and $\sigma$ $:{\mathbb R}^d\to{\mathbb R}^d\otimes{\mathbb R}^d$ satisfying {\bf (H$_\sigma$)} with some $p_0>d$}.
\end{align}
and consider the following SDE on probability space
$(\Omega,{\mathscr F},({\mathscr F}_t)_{t\ge0}, {\mathbb P})$:
\begin{equation}\label{SZZ001}
X_{s,t}^x=x+\int_s^tB(r,X_{s,r}^x){\mathord{{\rm d}}} r
+\int_s^t\sigma(X_{s,r}^x){\mathord{{\rm d}}} W_r,
\end{equation}
where $W_t$ is a standard $d$-dimensional Brownian motion and the PDE on $[0,T]\times{\mathbb R}^d$
\begin{align}\label{NBB00}
\partial_tu=a_{ij}\partial_i\partial_ju-\lambda u+B\cdot\nabla u+f,\quad u(0)=\varphi,
\end{align}
where $\lambda\ge0$, $f\in \widetilde{{\mathbb L}}^{p_0}(T)$, $\varphi\in C^\infty_b$ and $a_{ij}:=\frac{1}{2}\sum\limits_{k=1}^{d}\sigma_{ik}\sigma_{jk}$.
Under condition \eqref{NBB03}, by \cite{Xi-Xi-Zh-Zh}
there is a unique strong solution $X_{s,\cdot}^x$ to \eqref{SZZ001} for any $(s,x)\in{\mathbb R}_+\times{\mathbb R}^d$.
The purpose of this section is to obtain some moment estimates for the following functionals of $X_{0,t}^x$
\begin{align*}
\int_0^T f(s,X^x_{0,\pi_h(s)}){\mathord{{\rm d}}} s\quad \text{and}\quad\int_0^T\left[f(s,X^x_{0,s},\mu^x_{s})-f(s,X^x_{0,\pi_h(s)},\mu^x_{\pi_h(s)}) \right]{\mathord{{\rm d}}} s,
\end{align*}
where $\mu^x_t$ is the distribution of $X^x_{0,t}$ and $f\in {\mathbb L}^p_q(T)$ with some $2/q+d/p<2$.
The first one is classically based on Krylov's type estimates. Compared to the smooth coefficient case in \cite{SXW2022}, the $f$ in second one has no regularity. To overcome this obstacle, we use the observation mentioned in introduction to replace $f(X^x_{0,t})-f(X^x_{0,\pi_h(t)})$ by $P^X_{0,t}f-P^X_{0,\pi_h(s)}f$, where $P^X$ is the semi-group generated by $X$. Hence, we only need to obtain some time regularity results for the semigroup.

In subsection \ref{sec:3.1}, we consider the time-homogeneous case when $B\equiv0$. By Girsanov's Theorem, we extend the results in subsection \ref{sec:3.1} to $X^x_{s,t}$ in subsection \ref{Sec3.1}. Moreover, we obtain additional time regularity estimates for $P^X$ by Duhamel's formula which can not be gotten from Girsanov's Theorem. In the light of Duhamel's formula again, we also have two time regularity estimates for $\nabla u$ in subsection \ref{SubS32}, where $u$ is the solution to \eqref{NBB00}.

For simplicity, throughout this section we denote by
\begin{align*}
\Xi:=(d,T,p_0,|\mspace{-3mu}|\mspace{-3mu}| B|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{{\mathbb L}}}^{p_0}(T)},\kappa_1,\beta)
\end{align*}
the set of parameters.

}

\subsection{Time regularity for solutions to SDE with no drift}\label{sec:3.1}
First of all,
we introduce the following generalized It\^o's formula which is in \cite[Lemma 4.1-iii)]{Xi-Xi-Zh-Zh} (see also \cite[Theorem 3.7]{KR05} for an original one).
\begin{lemma}[Generalized It\^o's formula]
Let $p,q\in[2,\infty)$ with $2/q+d/p<1$.
For any $T>0$ and measurable function $u$ on ${\mathbb R}_+\times{\mathbb R}^d$
with $\partial_tu\in \widetilde{{\mathbb L}}^p_q(T)$ and $u\in\widetilde{{\mathbb H}}_q^{2,p}(T)$,
we have for any $t\in[s,T]$ and $x\in{\mathbb R}^d$,
\begin{align}\label{GIF001}
\begin{split}
u(t,X^x_{s,t})=u(s,x)&+\int_s^t(\partial_ru+a_{ij}\partial_i\partial_j u
+B\cdot\nabla u)(r,X_{s,r}(x)){\mathord{{\rm d}}} r\\
&+\int_s^t \nabla u(r,X_{s,r}(x)){\mathord{{\rm d}}} W_r.
\end{split}
\end{align}
\end{lemma}

In this subsection, we consider the following case when $B\equiv0$
\begin{equation}\label{SZZ00}
Z_t^x=x+\int_0^t\sigma(Z_s^x){\mathord{{\rm d}}} W_s.
\end{equation}
Define $P_t^\sigma f(x):={\mathbb E} f(Z_t^x)$. By Proposition \ref{well-PDE}, there is a unique solution to the following second order parabolic PDE in $[0,T]\times{\mathbb R}^d$:
\begin{align}\label{PreZ:PDE}
\partial_tu=a_{ij}\partial_i\partial_ju,\quad u_0=\varphi.
\end{align}

\begin{lemma}\label{PR-PDE:Z}
Assume {\bf (H$_\sigma$)} holds. Let $0\leqslant s\leqslant t$, $\varphi\in C_b^\infty(\R^d)$, $u$ and
$Z_{t}^x$ be the solution to \eqref{PreZ:PDE} and \eqref{SZZ00} respectively. Then we have ${\mathbb P}$-a.s.
\begin{equation}\label{PR-PDE:Zeq}
{\mathbb E}\left[\varphi(Z_t^x)|{\mathscr F}_s\right]=u(t-s,Z_s^x).
\end{equation}
In particular,
\begin{align}\label{ChH00}
{\mathbb E}\left[\varphi(Z_t^x)|{\mathscr F}_s\right]=P^\sigma_{t-s}\varphi(Z_s^x)\quad
{\mathbb P}-{\rm a.s.}
\end{align}
Moreover, for any $t\ge0$ and $f\in C^2_b$,
\begin{align}\label{NAA01}
P^\sigma_{t}f-f=\int_0^t P^\sigma_{r}(a_{ij}\partial_i\partial_jf){\mathord{{\rm d}}} r.
\end{align}
\end{lemma}
\begin{proof}
For all $t>0$, using the generalized It\^o's formula \eqref{GIF001}
to $s\mapsto u(t-s,Z_s^x)$, we have
\begin{align*}
u(0,Z_t^x)=u(t-s,Z_s^x)+\int_s^t\left(-\partial_ru(t-r,Z_r^x)
+a_{ij}\partial_i\partial_ju(t-r,Z_r^x)\right){\mathord{{\rm d}}} r
+\int_s^t\nabla u(t-r,Z_r^x){\mathord{{\rm d}}} W_r.
\end{align*}
Noting that $u(t-s,Z_s^x)$ is ${\mathscr F}_s$-measurable,
taking conditional expectations with respect to ${\mathscr F}_s$ on both sides, we have
\begin{align*}
{\mathbb E}\left[\varphi(Z_t^x)|{\mathscr F}_s\right]=u(t-s,Z_s^x)\quad
{\mathbb P}-{\rm a.s.},
\end{align*}
which implies that
\begin{align*}
P^\sigma_t\varphi(x)={\mathbb E}\varphi(Z_t^x)=u(t,x)
\end{align*}
provided by taking $s=0$. Then \eqref{ChH00} is straightforward from \eqref{PR-PDE:Zeq}.
For \eqref{NAA01}, since $f\in C^2_b$, we use the classical It\^o's formula and have
\begin{align*}
f(Z^x_{t})=f(x)+\int_0^t a_{ij}\partial_i\partial_j f(Z^x_{r}){\mathord{{\rm d}}} r+\int_0^t \nabla f(Z^x_{r}){\mathord{{\rm d}}} W_r.
\end{align*}
Then, we have \eqref{NAA01} by taking expectations and complete the proof.
\end{proof}

Based on Lemma \ref{PR-PDE:Z} and the uniqueness of \eqref{PreZ:PDE},
we have the following Chapman-Kolmogorov equation
\begin{align}\label{CKE00}
P^\sigma_{s}P^\sigma_{t}=P^\sigma_{s+t}.
\end{align}

Denote by
$$
\rho_t(x):=(2\pi t)^{-d/2}e^{-|x|^2/(2t)}
$$
the density of $W_t$. Then the following lemma is from \cite[Theorem 2.3]{CHXZ}.
\begin{lemma}\label{TDeqlem}
Assume {\bf (H$_\sigma$)} holds. There is a unique function
$p_\cdot^\sigma(\cdot,\cdot):\R_+\times\R^{2d}\rightarrow\R$ such that for any $j=0,1,2$
\begin{equation}\label{TDeq}
|\nabla_x^jp_t^\sigma(x,y)|\leqslant c_1t^{-\frac{j}{2}}\rho_{c_2t}(x-y)
\end{equation}
and
\begin{equation}\label{PR-RE}
  P_t^\sigma f(x)=\int_{\R^d}f(y)p_t^\sigma(x,y){\mathord{{\rm d}}} y
\end{equation}
for any $f\in C({\mathbb R}^d)$, where $c_1$ and $c_2>0$ are positive constants depending on $\Xi$.
\end{lemma}

For any $h\in(0,1)$, set $\pi_h(t):=t$ when $t\in[0,h)$ and
\begin{align*}
\pi_h(t):={kh},\ \ t\in[{kh},{(k+1)h}), ~k\ge1.
\end{align*}
\begin{remark}
The reason why we define $\pi_h(t)=t$ for $t\in[0,h)$ is that the function space here is $L^p$.
If the initial data don't have a $L^q$ density, ${\mathbb E} f(Z_{\pi_h(t)})={\mathbb E} f(Z_0)$, $f\in L^p$, will blow up for all $t<{h}$.
\end{remark}
Now we give the following Krylov's estimate and Khasminskii's estimate.
\begin{lemma}
Assume {\bf (H$_\sigma$)} holds. For any $T>0$, $k=0,1,2$, $p\in[1,\infty]$ and $q\in[p,\infty]$, there is a constant $C=C(\Xi,p,q)$ such that for all $0\leqslant s<t\leqslant T$, $x\in{\mathbb R}^d$ and positive function $f\in\widetilde{L}^p$
\begin{align}\label{BB00}
|\mspace{-3mu}|\mspace{-3mu}| \nabla^kP^\sigma_tf|\mspace{-3mu}|\mspace{-3mu}|_q\leqslant Ct^{-k/2-d/(2p)+d/(2q)}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}
\end{align}
and for $2/q+d/p<2$, $h>0$ and positive function $f\in\widetilde{{\mathbb L}}^p_q(T)$
\begin{align}\label{Kry00}
{\mathbb E}\int_s^t f(r,Z_r^x){\mathord{{\rm d}}} r+{\mathbb E}\int_s^t f(r,Z^x_{\pi_h(r)}){\mathord{{\rm d}}} r\leqslant C(t-s)^{1-\frac{1}{q}-\frac{d}{2p}}
|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}.
\end{align}
Moreover, for any $f\in\widetilde{{\mathbb L}}^p_q(T)$ with $d/p+2/q<2$,
\begin{align}\label{Kh00}
\sup_{x\in{\mathbb R}^d}{\mathbb E}\exp\left(\int_0^T f(t,Z^x_t){\mathord{{\rm d}}} t\right)<\infty.
\end{align}

\end{lemma}
\begin{proof}
Without loss of generality we assume that $c_2=1$ in \eqref{TDeq}.
Combining Lemma \ref{TDeqlem} and Young's convolution inequality, one sees that
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}|\nabla^kP^\sigma_{t}f|\mspace{-3mu}|\mspace{-3mu}|_{q}
&
\lesssim |\mspace{-3mu}|\mspace{-3mu}|\int_{\R^d}f(y)\nabla_x^kp_t^\sigma(\cdot,y){\mathord{{\rm d}}} y|\mspace{-3mu}|\mspace{-3mu}|_q\\
&
\lesssim t^{-k/2}|\mspace{-3mu}|\mspace{-3mu}|\rho_{t}*f|\mspace{-3mu}|\mspace{-3mu}|_q
\lesssim t^{-k/2}\sup_w\|{\mathbf{1}}_{|\cdot-w|\leqslant 1}\int_{{\mathbb R}^d}f(\cdot-y)\rho_t(y){\mathord{{\rm d}}} y\|_q\\
&
\lesssim t^{-k/2}\sup_w\|{\mathbf{1}}_{|\cdot-w|\leqslant 1}
\frac{1}{|B_1|}\int_{{\mathbb R}^d}\int_{{\mathbb R}^d}{\mathbf{1}}_{|y-z|\leqslant 1}
f(\cdot-y){\mathbf{1}}_{|y-z|\le1}\rho_t(y){\mathord{{\rm d}}} y{\mathord{{\rm d}}} z\|_q\\
&
\lesssim t^{-k/2} \sup_w\|\int_{{\mathbb R}^d}\int_{{\mathbb R}^d}{\mathbf{1}}_{|\cdot-y-w+z|\leqslant 2}f(\cdot-y){\mathbf{1}}_{|y-z|\le1}\rho_t(y){\mathord{{\rm d}}} y{\mathord{{\rm d}}} z\|_q\\
&
\lesssim t^{-k/2}\sup_w \int_{{\mathbb R}^d}\|\int_{{\mathbb R}^d}{\mathbf{1}}_{|\cdot-y-w+z|\leqslant 2}f(\cdot-y){\mathbf{1}}_{|y-z|\le1}\rho_t(y){\mathord{{\rm d}}} y\|_q{\mathord{{\rm d}}} z\\
&
\lesssim t^{-k/2}\int_{{\mathbb R}^d}\sup_w\|{\mathbf{1}}_{|\cdot-w+z|\leqslant 2}f(\cdot)\|_p\|{\mathbf{1}}_{|\cdot-z|\le1}\rho_t(\cdot)\|_r{\mathord{{\rm d}}} z\\
&
\lesssim  t^{-k/2} \int_{{\mathbb R}^d}\left(\int_{|y-z|\le1}|\rho_t(y)|^{r}{\mathord{{\rm d}}} y\right)^{1/r}{\mathord{{\rm d}}} z|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p,
\end{align*}
where $1+1/q=1/r+1/p$.
Next, one sees that
\begin{align*}
\int_{{\mathbb R}^d}\left(\int_{|y-z|\le1}(\rho_t(y))^{r}{\mathord{{\rm d}}} y\right)^{1/r}{\mathord{{\rm d}}} z
&
\lesssim \|\rho_t\|_{r}+\int_{|z|>2}\left(\int_{|y-z|\le1}(\rho_t(y))^{r}{\mathord{{\rm d}}} y\right)^{1/r}{\mathord{{\rm d}}} z.\\
\end{align*}
Noting that
\begin{align*}
|z|>2,|y-z|\le1\Rightarrow |y|\geqslant|z|-|y-z|\geqslant\frac{|z|}{2},
\end{align*}
which implies that $\rho_t(y)\lesssim \rho_{t}(z/2)$, and
$\|\rho_t\|_{r}\lesssim t^{-d/2+d/(2r)}=t^{-d/(2p)+d/(2q)}$, we have
\begin{align*}
\int_{{\mathbb R}^d}\left(\int_{|y-z|\le1}(\rho_t(y))^{r}{\mathord{{\rm d}}} y\right)^{1/r}{\mathord{{\rm d}}} z
&
\lesssim \|\rho_t\|_{r}+ \int_{|z|>2}\left(\int_{|y-z|\le1}(\rho_t(z/2))^{r}{\mathord{{\rm d}}} y\right)^{1/r}{\mathord{{\rm d}}} z\\
&
\lesssim t^{-d/(2p)+d/(2q)}+ \int_{|z|>2}\rho_t(z/2){\mathord{{\rm d}}} z\lesssim t^{-d/(2p)+d/(2q)}+1
\end{align*}
and obtain \eqref{BB00}.

Now we show \eqref{Kry00}. Set $1/p+1/p'=1/q+1/q'=1$. Without loss of generality, we take $s=0$. By \eqref{BB00}, for any $h>0$,
\begin{align*}
{\mathbb E}\int_0^t f(s,Z^x_{\pi_h(s)}){\mathord{{\rm d}}} s\lesssim \int_{{\mathbb R}^d}\Big(\int_0^t\(\int_{|y-z|\le1}|\rho_{\pi_h(s)}(y)|^{p'}{\mathord{{\rm d}}} y\)^{q'/p'}{\mathord{{\rm d}}} s\Big)^{1/q'}{\mathord{{\rm d}}} z|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}.
\end{align*}
Then, we have
\begin{align*}
{\mathscr I}&:=\int_{{\mathbb R}^d}\Big(\int_0^t\(\int_{|y-z|\le1}|\rho_{\pi_h(s)}(y)|^{p'}{\mathord{{\rm d}}} y\)^{q'/p'}{\mathord{{\rm d}}} s\Big)^{1/q'}{\mathord{{\rm d}}} z\\
&\lesssim \(\int_0^t\|\rho_{\pi_h(s)}\|_{p'}^{q'}{\mathord{{\rm d}}} s\)^{1/q'}+\int_{|z|>2}\Big(\int_0^t|\rho_{\pi_h(s)}(z/2)|^{q'}{\mathord{{\rm d}}} s\Big)^{1/q'}{\mathord{{\rm d}}} z\\
&\lesssim  \(\int_0^t(\pi_h(s))^{-dq'/(2p)}{\mathord{{\rm d}}} s\)^{1/q'}+t\int_{|z|>2}|z|^{-d}\exp(-|z|^2/(16T)){\mathord{{\rm d}}} z\\
&\lesssim t^{1-1/q-d/(2p)}+t
\end{align*}
provided by $dq'/(2p)<1$ and $\rho_s(z)\leqslant C|z|^{-d}\exp(-\frac{|z|^2}{4T})$ for all $s<T$.
Similarly, we obtain
$$
{\mathbb E}\int_0^tf(s,Z_s^x){\mathord{{\rm d}}} s\lesssim t^{1-\frac1q-\frac{d}{2p}}+t.
$$

Finally, noting that by \eqref{ChH00}
\begin{align*}
{\mathbb E}\Big[\int_s^t f(s,Z_r^x){\mathord{{\rm d}}} r\Big|{\mathscr F}_s\Big]={\mathbb E}\int_s^t f(s,Z_{r-s}^y){\mathord{{\rm d}}} r\Big|_{y=Z_s^x},
\end{align*}
 \eqref{Kh00} is direct from \eqref{Kry00} (see \cite[Corollary 3.5]{Zh-Zh1} for example) and we complete the proof.

\end{proof}

\begin{remark}
Notice that for any fixed $h>0$ and $x\in{\mathbb R}^d$
\begin{align*}
{\mathbb E}\exp\(\int_0^T f(t,Z^x_{\pi_h(t)}){\mathord{{\rm d}}} t\)<\infty
\end{align*}
is not true. For example, letting $Z^x_t=W_t$, when $d\ge2$ and $f(t,x)=|x|^{-1/2}\in \tilde{L}^{d+1}({\mathbb R}^d)$, we have
\begin{align*}
&\quad{\mathbb E}\exp\left(\int_{h}^{2h} f(t,W_{\pi_h(t)}){\mathord{{\rm d}}} t\right)={\mathbb E}\exp(h|W_{h}|^{-1/2})\\
&=\int_{{\mathbb R}^d}e^{\frac{h}{\sqrt{|x|}}}\rho_h(x){\mathord{{\rm d}}} x\geqslant\frac{h^{2d}}{(2d)!}\int_{{\mathbb R}^d}\frac{1}{|x|^d}\rho_h(x){\mathord{{\rm d}}} x=\infty.
\end{align*}

\end{remark}

\begin{lemma}\label{lemmm03}
Let $T>0$, $k=0,1$, and $1\leqslant p\leqslant q\leqslant\infty$. There is a constant $C=C(\Xi,p,q)$ such that for all $0< s\leqslant t \leqslant T$,
\begin{equation}\label{MM03}
|\mspace{-3mu}|\mspace{-3mu}| \nabla^k(P^\sigma_t\varphi-P^\sigma_s\varphi)|\mspace{-3mu}|\mspace{-3mu}|_q\leqslant C\Big([(t-s)^{\frac{2-k}{2}}s^{\frac{k-2}{2}}]\wedge1\Big)s^{-\frac{k}{2}-\frac{d}{2p}
+\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_p.
\end{equation}

\end{lemma}

\begin{proof}
Based on \eqref{CKE00} and \eqref{NAA01}, one sees that
\begin{align*}
P^\sigma_t\varphi-P^\sigma_s\varphi=P^\sigma_{t-s}(P^\sigma_s\varphi)-P^\sigma_s\varphi=\int_0^{t-s}P^\sigma_r(a_{ij}\partial_i\partial_jP^\sigma_s\varphi){\mathord{{\rm d}}} r.
\end{align*}
By \eqref{BB00}, we have
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| \nabla^k(P^\sigma_t\varphi-P^\sigma_s\varphi)|\mspace{-3mu}|\mspace{-3mu}|_q
&
\lesssim\int_0^{t-s}r^{-\frac{k}{2}}|\mspace{-3mu}|\mspace{-3mu}|\nabla^2P_s^\sigma\varphi|\mspace{-3mu}|\mspace{-3mu}|_q{\mathord{{\rm d}}} r\\
&
\lesssim\int_0^{t-s}r^{-\frac{k}{2}}s^{-1+\frac{d}{2q}-\frac{d}{2p}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_p{\mathord{{\rm d}}} r\\
&
\lesssim \Big[(t-s)^{\frac{2-k}{2}}s^{\frac{k-2}{2}}\Big]s^{-\frac{k}{2}-\frac{d}{2p}
+\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_p.
\end{align*}
Moreover, noting that
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| \nabla^k(P^\sigma_t\varphi-P^\sigma_s\varphi)|\mspace{-3mu}|\mspace{-3mu}|_q&\leqslant |\mspace{-3mu}|\mspace{-3mu}| \nabla^kP^\sigma_t\varphi|\mspace{-3mu}|\mspace{-3mu}|_q+|\mspace{-3mu}|\mspace{-3mu}|\nabla^k P^\sigma_s\varphi|\mspace{-3mu}|\mspace{-3mu}|_q\\
&\lesssim s^{-\frac{k}{2}-\frac{d}{2p}
+\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_p
\end{align*}
provided by \eqref{BB00} and $s\leqslant t$, we complete the proof.
\end{proof}

Now, we introduce the following results for $Z^x_t$. When $p=\infty$ and $\sigma\equiv{\mathbb I}$, the following lemma has been proved in \cite[Lemma 2.1]{DG20} for Brownian motion. For classical $L^p$ space, L\^e and Ling obtained these results by stochastic sewing lemma in \cite{LL21}. For the localized $L^p$ space, we provide a different proof here,
which is based on \eqref{MM03} and \eqref{BB00}.
\begin{lemma}\label{CruL}
Assume {\bf (H$_\sigma$)} holds. For any $T>0$ and $p\in(d\vee 2,\infty)$, there is a constant $C=C(\Xi,p)$ such that for any stopping time $\tau\leqslant T$,  $h\in(0, 1)$, $x\in{\mathbb R}^d$ and $f\in\widetilde{{\mathbb L}}^p(T)$,
\begin{align}\label{BB03}
{\mathbb E}\Big|\int_{0}^\tau(f(t,Z_t^x)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\Big|^2
\leqslant C h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2.
\end{align}

\end{lemma}
\begin{proof}
We divided the proof into four steps. In step 1, we prove that
\begin{align}\label{MM00}
{\mathbb E}\Big|\int_{s}^t(f(r,Z^x_r)-f(r,Z^x_{\pi_h(r)})){\mathord{{\rm d}}} r\Big|^2
\leqslant C \(h s^{-\frac{d}{2p}}(t-s)^{1-\frac{d}{2p}}+ h\ln h^{-1}(t-s)^{1-\frac{d}{p}}\)|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2;
\end{align}
 In step 2, we show \eqref{BB03} for $\tau=T$; In step 3, we show
\begin{align}\label{DD000}
\sup_{a,b\in[0,T]}{\mathbb E}\left(\int_{a}^b(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right)^2
\leqslant C h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2;
\end{align}
In step 4, we show \eqref{BB03} for any stopping time $\tau$ with $\tau\leqslant T$.\\
\vspace{2mm}

{\bf (Step 1)}
Without loss of generality, we drop the time $t$ in $f(t,x)$. Firstly, notice that by H\"older's inequality and \eqref{Kry00},
\begin{align*}
&\quad{\mathbb E}\left(\int_0^{{2h}}f(Z^x_t)-f(Z^x_{\pi_h(t)}){\mathord{{\rm d}}} t\right)^2\le2h{\mathbb E}\int_0^{{2h}}\Big|f(Z^x_t)-f(Z^x_{\pi_h(t)})\Big|^2{\mathord{{\rm d}}} t\lesssim h|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p.
\end{align*}
Hence, without loss of generality, we assume $s>2h$.
Based on the symmetric,
\begin{align*}
&
\quad{\mathbb E}\Big|\int_{s}^t(f(Z^x_r)-f(Z^x_{\pi_h(r)})){\mathord{{\rm d}}} r\Big|^2\\
=&2\int_{s}^t\int_{{r_1}}^t {\mathbb E}\[(f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)}))(f(Z^x_{r_2})-f(Z^x_{\pi_h(r_2)}))\]{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
=&2\int_{s}^t\int_{{r_1}}^{r_1+h} {\mathbb E}\[(f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)}))(f(Z^x_{r_2})-f(Z^x_{\pi_h(r_2)}))\]{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&+2\int_{s}^t\int_{{r_1+h}}^t {\mathbb E}\[(f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)}))(f(Z^x_{r_2})-f(Z^x_{\pi_h(r_2)}))\]{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
:=&2{\mathscr I}_1+2{\mathscr I}_2.
\end{align*}
By H\"older's inequality and \eqref{BB00}, one sees that
\begin{align*}
&
\quad{\mathbb E}\[(f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)}))(f(Z^x_{r_2})-f(Z^x_{\pi_h(r_2)}))\]\\
&
\leqslant \Big({\mathbb E}\Big|f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)})\Big|^2\Big)^{1/2}\Big({\mathbb E}\Big|f(Z^x_{r_2})-f(Z^x_{\pi_h(r_2)})\Big|^2\Big)^{1/2}\\
&
\lesssim (r_1^{-\frac{d}{2p}}+{\pi_h(r_1)}^{-\frac{d}{2p}})(r_2^{-\frac{d}{2p}}+{\pi_h(r_2)}^{-\frac{d}{2p}})|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p^2,
\end{align*}
which implies that
\begin{align*}
{\mathscr I}_1&\lesssim|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p^2 \int_{s}^t\int_{{r_1}}^{r_1+h}(r_1^{-\frac{d}{2p}}+{\pi_h(r_1)}^{-\frac{d}{2p}})(r_2^{-\frac{d}{2p}}+{\pi_h(r_2)}^{-\frac{d}{2p}}){\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p^2(s-h)^{-\frac{d}{2p}} h \int_{s}^t(r_1^{-\frac{d}{2p}}+{\pi_h(r_1)}^{-\frac{d}{2p}}){\mathord{{\rm d}}} r_1\lesssim hs^{-\frac{d}{2p}} (t-s)^{1-\frac{d}{2p}}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p^2.
\end{align*}

For ${\mathscr I}_2$, we use conditional expectation and Markov property \eqref{ChH00}. For simplicity, let
\begin{align*}
{\mathbb E}^{{\mathscr G}}[\cdot]:={\mathbb E}[\cdot|{\mathscr G}].
\end{align*}
Then, noting that $r_1\leqslant r_2-h<r_2$, by Markov property \eqref{ChH00}, we have
\begin{align*}
{\mathscr I}_2&=\int_{s}^t\int_{{r_1+h}}^t {\mathbb E}\Big((f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)})){\mathbb E}^{{\mathscr F}_{r_1}}[f(Z^x_{r_2})-f(Z^x_{\pi_h(r_2)})]\Big){\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&=\int_{s}^t\int_{{r_1+h}}^t {\mathbb E}\Big(\[f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)})\]\[P^\sigma_{r_2-r_1}f(Z^x_{r_1})-P^\sigma_{\pi_h(r_2)-r_1}f(Z^x_{r_1})\]\Big){\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1.
\end{align*}
Setting $G:=P^\sigma_{r_2-r_1}f-P^\sigma_{\pi_h(r_2)-r_1}f$, in view of H\"older's inequality and \eqref{BB00}, one sees that
\begin{align*}
{\mathscr I}_2&\leqslant \int_{s}^t\int_{{r_1+h}}^t \Big({\mathbb E}\Big|f(Z^x_{r_1})-f(Z^x_{\pi_h(r_1)})\Big|^2\Big)^{1/2}\Big({\mathbb E} |G(Z^x_{r_1})|^2\Big)^{1/2}{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&\lesssim |\mspace{-3mu}|\mspace{-3mu}| |f|^2|\mspace{-3mu}|\mspace{-3mu}|_{p/2}^{1/2}\int_{s}^t\int_{{r_1+h}}^t ({r_1}^{-\frac{d}{2p}}+{\pi_h(r_1)}^{-\frac{d}{2p}}){r_1}^{-\frac{d}{2p}}|\mspace{-3mu}|\mspace{-3mu}| |G|^2|\mspace{-3mu}|\mspace{-3mu}|_{p/2}^{1/2}{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}\int_{s}^t\int_{{r_1+h}}^t {(r_1-h)}^{-\frac{d}{p}}|\mspace{-3mu}|\mspace{-3mu}| G|\mspace{-3mu}|\mspace{-3mu}|_{p}{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1.
\end{align*}
Noting that by \eqref{MM03},
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| G|\mspace{-3mu}|\mspace{-3mu}|_{p}&=|\mspace{-3mu}|\mspace{-3mu}| P^\sigma_{r_2-r_1}f-P^\sigma_{\pi_h(r_2)-r_1}f|\mspace{-3mu}|\mspace{-3mu}|_p\\
&\lesssim \Big([(r_2-\pi_h(r_2))(\pi_h(r_2)-r_1)^{-1}]\wedge 1\Big)|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p\\
&\lesssim \Big[(h(r_2-r_1-h)^{-1})\wedge 1\Big]|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p,
\end{align*}
 we have
\begin{align*}
{\mathscr I}_2
&
\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}^2\int_{s}^t\int_{{r_1+h}}^t {(r_1-h)}^{-\frac{d}{p}}\Big[(h(r_2-r_1-h)^{-1})\wedge 1\Big]{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&
\lesssim|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}^2\int_{s-h}^{t-h}\int_{0}^{t} {(r_1)}^{-\frac{d}{p}}\Big[(h (r_2)^{-1})\wedge 1\Big]{\mathord{{\rm d}}} r_2{\mathord{{\rm d}}} r_1\\
&
\lesssim h(t-s)^{1-\frac{d}{p}}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}^2\int_{0}^{t/h}\Big[(r_2)^{-1}\wedge 1\Big]{\mathord{{\rm d}}} r_2\lesssim h\ln h^{-1}(t-s)^{1-\frac{d}{p}} |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}^2,
\end{align*}
provided by a change of variable, and we obtain \eqref{MM00}.\\

{\bf (Step 2)} In this step, we use the method in \cite[Lemma 3.4]{LL21} to show \eqref{BB03} with $\tau=T$. Let $\alpha:=\frac{d}{4p}$, $\|\cdot\|:=\({\mathbb E}|\cdot|^2\)^{\frac12}$,
\begin{align*}
F_t:=\int_{0}^t(f(r,Z^x_r)-f(r,Z^x_{\pi_h(r)})){\mathord{{\rm d}}} r\quad \text{and} \quad {\mathcal H}:=\(C h\ln h^{-1}\)^{1/2}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)},
\end{align*}
where $C$ is the constant in \eqref{MM00}.
By \eqref{MM00}, we have
\begin{align}\label{MM01}
\|F_t-F_s\|\leqslant {\mathcal H} (s^{-\alpha}(t-s)^{1/2-\alpha}+(t-s)^{1/2-2\alpha}).
\end{align}
For any $n\in{\mathbb N}$, set
\begin{align*}
t_n:=2^{-n}T.
\end{align*}
Then, by \eqref{MM01}
\begin{align*}
\|F_T\|\leqslant\sum_{n=0}^\infty\|F_{t_n}-F_{t_{n+1}}\|\leqslant{\mathcal H}\sum_{n=0}^\infty (t_{n+1}^{-\alpha}(t_{n}-t_{n+1})^{1/2-\alpha}+(t_n-t_{n+1})^{1/2-2\alpha}).
\end{align*}
Noting that $t_{n+1}=T2^{-n-1}$ and $t_n-t_{n+1}=T2^{-n-1}$, we have
\begin{align*}
\|F_T\|
&
\leqslant{\mathcal H}\sum_{n=0}^\infty (T^{-\alpha}2^{\alpha(n+1)}(t_{n}-t_{n+1})^{1/2-\alpha}
+(t_n-t_{n+1})^{1/2-2\alpha})\\
&
\leqslant 2{\mathcal H} T^{1/2-2\alpha}\sum_{n=0}^\infty 2^{-(1/2-2\alpha)(n+1)}\leqslant 2{\mathcal H} T^{1/2-2\alpha}
\end{align*}
and obtain \eqref{BB03} for $\tau=T$.

{\bf (Step 3)}
For any $a,b\in[0,T]$, define
\begin{align*}
f_{a,b}(t,x):={\mathbf{1}}_{t\in[a,b]}f(t,x).
\end{align*}
 Then, by Step 1, one sees that
 \begin{align*}
{\mathbb E}\left(\int_{a}^b(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right)^2
&
={\mathbb E}\(\int_{0}^T(f_{a,b}(t,Z^x_t)-f_{a,b}(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\)^2\\
&
\leqslant C h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f_{a,b}|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2.
\end{align*}
Noting that $|\mspace{-3mu}|\mspace{-3mu}| f_{a,b}|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}\leqslant |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}$, we get \eqref{DD000}.

{\bf (Step 4)}
Without loss of generality, we assume that $h<T/2$ and $\tau$ only takes finite value $a_1,a_2,...,a_n\in[0,T]$. Actually, for any stopping time, we choose $\tau_n$ which takes finite value to approximate $\tau$ and \eqref{BB03} comes from \eqref{Kry00} and dominated convergence theorem. Firstly,
\begin{align*}
{\mathbb E}\left|\int_{\tau}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right|^2
&
\leqslant 2{\mathbb E}\left|\int_{\tau}^{(\tau+2h)\wedge T}(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right|^2\\
&\quad
+2{\mathbb E}\left|\int_{(\tau+2h)\wedge T}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right|^2.
\end{align*}
Notice that by H\"older's inequality and \eqref{Kry00} for $p/2>d/2$,
\begin{align*}
&
\quad{\mathbb E}\Big|\int_{\tau}^{(\tau+2h)\wedge T}(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\Big|^2\\
&
\leqslant {\mathbb E}\Big|\int_{0}^{T}{\mathbf{1}}_{t\in[\tau,\tau+2h]}{\mathord{{\rm d}}} t\int_{0}^{T}|f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})|^2{\mathord{{\rm d}}} t\Big|\\
&
\lesssim h |\mspace{-3mu}|\mspace{-3mu}| |f|^2|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p/2}(T)}\lesssim h |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p}(T)}^2.
\end{align*}
Now we estimate the second term. In fact,
\begin{align*}
&\quad
{\mathbb E}\left(\int_{(\tau+2h)\wedge T}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right)^2\\
&
=\sum_{i=1}^n{\mathbb E}\left[{\mathbf{1}}_{\tau=a_i}\left(\int_{(a_i+2h)\wedge T}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right)^2\right].
\end{align*}
Noting that ${\mathbf{1}}_{\tau=a_i}\in{\mathscr F}_{a_i}\subset{\mathscr F}_{([a_i/h]+1)h}$, without loss of generality, assuming $a_i+2h<T$, one sees that
\begin{align*}
&\quad
{\mathbb E}\left[{\mathbf{1}}_{\tau=a_i}\left(\int_{a_i+2h}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right)^2\right]\\
&
={\mathbb E}\Big({\mathbf{1}}_{\tau=a_i}{\mathbb E}\Big[\(\int_{a_i+2h}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\)^2\Big|{\mathscr F}_{a_i}\Big]\Big)\\
&
:={\mathbb E}({\mathbf{1}}_{\tau=a_i}{\mathscr A}_i),
\end{align*}
where
\begin{align*}
{\mathscr A}_i={\mathbb E}\left[{\mathbb E}\left[\left(\int_{a_i+2h}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right)^2\Big|{\mathscr F}_{([a_i/h]+1)h}\right]\Big|{\mathscr F}_{a_i}\right].
\end{align*}
Moreover, by Markov property \eqref{ChH00}, we have
\begin{align*}
&\quad
{\mathbb E}\Big[\(\int_{a_i+2h}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\)^2\Big|{\mathscr F}_{([a_i/h]+1)h}\Big]\\
&
={\mathbb E}\Big(\int_{a_i+2h}^T(f(t,Z^y_{t-([a_i/h]+1)h})-f(t,Z^y_{\pi_h(t)-([a_i/h]+1)h})){\mathord{{\rm d}}} t\Big)^2\Big|_{y=Z^x_{_{([a_i/h]+1)h}}}\\
&
={\mathbb E}\Big(
\int_{a_i+2h-([a_i/h]+1)h}^{T-([a_i/h]+1)h}
(f(t+([a_i/h]+1)h,Z^y_{t})-f(t+([a_i/h]+1)h,Z^y_{\pi_h(t)})){\mathord{{\rm d}}} t\Big)^2\Big|_{y=Z^x_{_{([a_i/h]+1)h}}}\\
&
\lesssim h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2,
\end{align*}
provided by \eqref{DD000}. Therefore, we have
\begin{align*}
{\mathbb E}\Big|\int_{\tau}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\Big|^2
&
\lesssim h |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2+\sum_{i=1}^n{\mathbb E}({\mathbf{1}}_{\tau=a_i}{\mathscr A}_i)\\
&
\lesssim  h |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2+h\ln h^{-1}
|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2\sum_{i=1}^n{\mathbb E}{\mathbf{1}}_{\tau=a_i}\\
&
\lesssim h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2,
\end{align*}
which implies that
\begin{align*}
&\quad
{\mathbb E}\left|\int_0^{\tau}(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right|^2\\
&
\leqslant 2{\mathbb E}\left|\int_{0}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right|^2
+2{\mathbb E}\left|\int_{\tau}^T(f(t,Z^x_t)-f(t,Z^x_{\pi_h(t)})){\mathord{{\rm d}}} t\right|^2\\
&
\lesssim h \ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2
\end{align*}
and complete the proof.

\end{proof}

\begin{corollary}
Assume {\bf (H$_\sigma$)} holds. For any $T>0$, $p\in(d\vee 2,\infty)$ and $\delta>0$, there is a constant $C=C(\Xi,p)$ such that for any $x\in{\mathbb R}^d$ and $f\in\widetilde{{\mathbb L}}^p(T)$,
\begin{align}\label{DD03}
{\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_{0}^t(f(s,Z^x_s)-f(s,Z^x_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)
\leqslant C h^{1-\delta}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2.
\end{align}

\end{corollary}
\begin{proof}
Let
\begin{align*}
\eta_t:=\Big|\int_{0}^t(f(s,Z^x_s)-f(s,Z^x_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2
\end{align*}
and $\eta^*_t:=\sup\limits_{s\in[0,t]}\eta_s$.
First of all, it follows from H\"older's inequality and \eqref{Kry00} that for any $\gamma\in(1,p/(d\vee 2))$,
\begin{align}
{\mathbb E}(\eta^*_T)^{\gamma}
&
\lesssim {\mathbb E}\int_{0}^T\left(|f|^{2\gamma}(s,Z^x_s)+|f|^{2\gamma}(s,Z^x_{\pi_h(s)})\right){\mathord{{\rm d}}} s\nonumber\\
&
\lesssim |\mspace{-3mu}|\mspace{-3mu}| |f|^{2\gamma}|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p/(2\gamma)}(T)}
\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^{2\gamma}.\label{DD00}
\end{align}
For any $\lambda>0$, let
\begin{align*}
\tau_\lambda:=\inf\{t\ge0,~\eta_t>\lambda\}.
\end{align*}
Notice that $\eta_{\tau_\lambda}=\lambda$ since $\eta$ is a continuous process. Then,
\begin{align*}
\lambda{\mathbb P}(\eta^*_T>\lambda)\leqslant\lambda{\mathbb P}(\tau_\lambda\leqslant T)
\leqslant {\mathbb E}\(\eta_{\tau_\lambda}{\mathbf{1}}_{\{\tau_\lambda\leqslant T\}}\)\leqslant {\mathbb E}\eta_{\tau_\lambda\wedge T}.
\end{align*}
In view of \eqref{BB03}, we have
\begin{align*}
\lambda{\mathbb P}(\eta^*_T>\lambda)\lesssim h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2.
\end{align*}
Set $\Xi_h:=h\ln h^{-1}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2$.
Hence, for any $\delta\in(0,1)$, by a change of variable,
\begin{align}
{\mathbb E}\(\eta^*_T\)^{1-\delta}
&
=(1-\delta)\int_0^\infty\lambda^{-\delta}
{\mathbb P}(\eta^*_T>\lambda){\mathord{{\rm d}}}\lambda\nonumber\\
&
\lesssim \int_0^\infty\lambda^{-\delta}
\(1\wedge(\Xi_h\lambda^{-1})\){\mathord{{\rm d}}}\lambda\nonumber\\
&
\lesssim \Xi_h^{1-\delta}\int_0^\infty \lambda^{-\delta}(1\wedge \lambda^{-1}){\mathord{{\rm d}}}\lambda
\lesssim (h\ln h^{-1})^{(1-\delta)}
\|f\|_{\widetilde{{\mathbb L}}^p(T)}^{2(1-\delta)}.\label{DD00000}
\end{align}
Combining \eqref{DD00} and \eqref{DD00000}, in view of H\"older's inequality, for any $\delta>0$ small enough, we have
\begin{align*}
{\mathbb E}\eta^*_T
&
={\mathbb E}\left[\(\eta^*_T\)^{1-\delta-\sqrt{\delta}}\(\eta^*_T\)^{\delta+\sqrt{\delta}}\right]\\
&
\lesssim \left[{\mathbb E}\(\eta^*_T\)^{1-\delta}\right]^{\frac{1-\delta-\sqrt{\delta}}{1-\delta}}
\Big[{\mathbb E}\(\eta^*_T\)^{(\sqrt{\delta}+1)(1-\delta)}\Big]^{\frac{\sqrt{\delta}}{1-\delta}}\\
&
\lesssim h^{1-\delta-\sqrt{\delta}}(\ln h^{-1})^{1-\delta-\sqrt{\delta}}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^{2}
\end{align*}
and complete the proof.

\end{proof}

\subsection{Time discretization for SDEs with $\widetilde{{\mathbb L}}^p$ drift}\label{Sec3.1}
Now, let's extend estimate \eqref{DD03} from $Z^x_{t-s}$ to the solution $X^x_{s,t}$ of SDE \eqref{SZZ001} both in the sense of path and distribution (see \eqref{BB05} and \eqref{LL02} below).
Recall that
\begin{align}\label{II03}
X^x_{s,t}=x+\int_s^t B(r,X^x_{s,r}){\mathord{{\rm d}}} r+\int_s^t\sigma(X^x_{s,r}){\mathord{{\rm d}}} W_r.
\end{align}
Denote by $\mu_{s,t}^x$ the distribution of $X^x_{s,t}$. For simplicity, we also denote by
\begin{align*}
X^x_t:=X^x_{0,t}\quad \text{and}\quad\mu_t^x:=\mu_{0,t}^x.
\end{align*}
The following estimates are from Girsanov's transform and estimates for $Z^x_{t}$ (see \cite[Lemma 4.1]{Xi-Xi-Zh-Zh} for (ii) and (iii)).
\begin{lemma}\label{Lem34}
Assume \eqref{NBB03} with $p_0>d\vee 2$.\\
(i) For any $T>0$ and $p\in(1,\infty)$, there is a constant $C=C(\Xi,p)$
such that for any $x\in{\mathbb R}^d$, $0\leqslant s<t\leqslant T$ and positive $f\in \widetilde{L}^p$,
\begin{align}\label{BB06}
{\mathbb E} f(X^x_{t})\leqslant C t^{-d/(2p)}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{p}.
\end{align}
(ii)For any $T>0$ and $2/q+d/p<2$, there is a constant $C=C(\Xi,p,q)$ such that for any $x\in{\mathbb R}^d$, $0\leqslant s<t\leqslant T$ and positive $f\in \widetilde{{\mathbb L}}^p_q(T)$,
\begin{align}\label{Kry02}
{\mathbb E}\int_s^t f(r,X^x_{r}){\mathord{{\rm d}}} r+{\mathbb E}\int_{s}^t f(r,X^x_{\pi_h(r)}){\mathord{{\rm d}}} r\leqslant C(t-s)^{1-\frac{1}{q}-\frac{d}{2p}}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p_q(T)}.
\end{align}
(iii)For any $T>0$, $d/p+2/q<2$ and $f\in \widetilde{{\mathbb L}}^p_q(T)$,
\begin{align}\label{Kh02}
\sup_{x\in{\mathbb R}^d}{\mathbb E}\exp\left(\int_0^T f(t,X^x_t){\mathord{{\rm d}}} t\right)<\infty.
\end{align}
(iv)For any $T>0$, $\delta>0$ and $p>d\vee2$, there is a constant $C=C(\Xi,p)$ such that for any $x\in{\mathbb R}^d$, $h\in(0,1)$ and $f\in \widetilde{{\mathbb L}}^p(T)$,
\begin{align}\label{BB05}
{\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_{0}^t(f(s,X^x_{s})-f(s,X^x_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)
\leqslant C h^{1-\delta}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^p(T)}^2.
\end{align}
\end{lemma}
\begin{proof}
Let $\widetilde{Z}^x$ be a solution on probability $(\widetilde{\Omega},\widetilde{{\mathscr F}},(\widetilde{{\mathscr F}}_t)_{t\ge0},\widetilde{{\mathbb P}})$ to the following SDE
\begin{align*}
\widetilde{Z}^x_t=x+\int_0^t\sigma(\widetilde{Z}^x_r){\mathord{{\rm d}}} \widetilde{W}_r,
\end{align*}
where $\widetilde{W}_t$ is a standard $d$-dimensional Brownian motion.
Since $p_0>2\vee d$, $B^2\in \widetilde{{\mathbb L}}^{p_0/2}(T)$. By \eqref{Kh00}, one sees that for any $\gamma>0$
\begin{align*}
\sup_x\widetilde{{\mathbb E}}\exp\left(\gamma\int_0^T |\sigma^{-1}B(t,\widetilde{Z}^x_t)|^2{\mathord{{\rm d}}} t\right)<\infty.
\end{align*}
Hence, by Novikov criterion,
\begin{align*}
{\mathbb Z}_T:= \exp\left(-\int_0^T \sigma^{-1}B(t,\widetilde{Z}^x_t){\mathord{{\rm d}}} \widetilde{W}_t-\frac{1}{2}\int_0^T|\sigma^{-1}B(t,\widetilde{Z}^x_t)|^2{\mathord{{\rm d}}} t\right)
\end{align*}
is integrable and for any $q>0$
\begin{align}\label{BB061}
\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}|{\mathbb Z}_T|^q\leqslant C(\Xi,q).
\end{align}
Define ${\mathord{{\rm d}}}{\mathbb Q}:={\mathbb Z}_T{\mathord{{\rm d}}} \widetilde{{\mathbb P}}$. Then, by Girsanov's theorem,
\begin{align*}
\text{$\bar{W}_t:=\widetilde{W}_t-\int_0^t \sigma^{-1}B(s,\widetilde{Z}^x_s){\mathord{{\rm d}}} s$ is a ${\mathbb Q}$-martingale}.
\end{align*}
In other words,
\begin{align*}
\widetilde{Z}^x_t=x+\int_0^t B(s,\widetilde{Z}_s){\mathord{{\rm d}}} s+\int_0^t\sigma(\widetilde{Z}_s){\mathord{{\rm d}}} \bar{W}_s,\ \ {\mathbb Q}-a.e.
\end{align*}
Therefore, by the uniqueness of \eqref{II03}, we have
\begin{align}\label{Tra001}
{\mathbb Q}\circ(\widetilde{Z}^x_\cdot)^{-1}={\mathbb P}\circ(X^x_\cdot)^{-1}.
\end{align}
Now, we show (i)-(iv) one by one.

(i) In view of \eqref{Tra001}, one sees that
\begin{align*}
{\mathbb E} f(X_t^x)=\widetilde{{\mathbb E}}_{{\mathbb Q}}f(\widetilde{Z}^x_t)=\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}\[{\mathbb Z}_Tf(\widetilde{Z}_t)\].
\end{align*}
By H\"older's inequality, \eqref{BB061} and \eqref{BB00}, we have for any $r\in(1,p)$ and $1/r'+1/r=1$,
\begin{align*}
{\mathbb E} f(X^x_t)\leqslant \left(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}|{\mathbb Z}_T|^{r'}\right)^{1/r'}\left(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}} |f(\widetilde{Z}_t)|^r\right)^{1/r}
\lesssim \left(t^{-dr/(2p)}|\mspace{-3mu}|\mspace{-3mu}| |f|^r|\mspace{-3mu}|\mspace{-3mu}|_{p/r}\right)^{1/r}
\lesssim t^{-d/(2p)}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_p,
\end{align*}
which is \eqref{BB06}.

(ii) Similar to (i), by H\"older's inequality and \eqref{Kry00}, we have
\begin{align*}
{\mathbb E} \int_s^tf(u,X^x_u){\mathord{{\rm d}}} u
&
\leqslant \left(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}|{\mathbb Z}_T|^{r'}\right)^{1/r'}\left(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}} \left[\int_s^tf(u,\widetilde{Z}^x_u){\mathord{{\rm d}}} u\right]^r\right)^{1/r}\\
&
\lesssim\left(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}\left( (t-s)^{r-1}\int_s^t|f(u,\widetilde{Z}^x_u)|^r{\mathord{{\rm d}}} u\right)\right)^{1/r}\\
&
\lesssim (t-s)^{1-1/q-d/(2p)}\left(|\mspace{-3mu}|\mspace{-3mu}| |f|^r|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p/r}_{q/r}}\right)^{1/r}\lesssim (t-s)^{1-1/q-d/(2p)}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p}_{q}}.
\end{align*}
The term ${\mathbb E} \int_s^tf(u,X^x_{\pi_h(u)}){\mathord{{\rm d}}} u$ is totally the same.

(iii)
For \eqref{Kh02}, it follows from H\"older's inequality again that
\begin{align*}
{\mathbb E}\exp\left(\int_0^T f(t,X^x_t){\mathord{{\rm d}}} t\right)
&
=\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}\left[{\mathbb Z}_T\exp\(\int_0^T f(t,\widetilde{Z}^x_t){\mathord{{\rm d}}} t\)\right]\\
&
\leqslant \(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}|{\mathbb Z}_T|^{r'}\)^{1/r'}\left(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}
\exp\(r\int_0^T f(t,\widetilde{Z}^x_t){\mathord{{\rm d}}} t\)\right)^{1/r}<\infty,
\end{align*}
provided by \eqref{BB061} and \eqref{Kh00}.

(iv) Let
\begin{align*}
A^f(h,X):=\sup_{t\in[0,T]}\Big|\int_{0}^t
(f(s,X^x_s)-f(s,X^x_{\pi_h(s)})){\mathord{{\rm d}}} u\Big|^2
\end{align*}
and
\begin{align*}
A^f(h,\widetilde{Z}):=\sup_{t\in[0,T]}\Big|\int_{{0}}^t
(f(s,\widetilde{Z}^x_s)-f(s,\widetilde{Z}^x_{\pi_h(s)})){\mathord{{\rm d}}} u\Big|^2.
\end{align*}
For any $\delta\in(0,1)$, notice that
\begin{align*}
{\mathbb E} A^f(h,X)
&
=\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}\({\mathbb Z}_TA^f(h,\widetilde{Z})\)\\
&
=\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}\({\mathbb Z}_T|A^f(h,\widetilde{Z})|^\delta|A^f(h,\widetilde{Z})|^{1-\delta}\).
\end{align*}
Based on H\"older's inequality, \eqref{BB061}, \eqref{Kry02} and \eqref{DD03}, for $1/r'+1/r=1$ with some $r\in(1,p/2)$, we have
\begin{align*}
{\mathbb E} A^f(h,X)
&
\leqslant \[\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}\(|{\mathbb Z}_T|^{1/\delta}A^f(h,\widetilde{Z})\)\]^{\delta}
\[\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}A^f(h,\widetilde{Z})\]^{1-\delta}\\
&
\leqslant\(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}|{\mathbb Z}_T|^{r'/\delta}\)^{\delta/r'}
\(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}|A^f(h,\widetilde{Z})|^{r}\)^{\delta/r}
\(\widetilde{{\mathbb E}_{\widetilde{{\mathbb P}}}}A^f(h,\widetilde{Z})\)^{1-\delta}\\
&
\lesssim |\mspace{-3mu}|\mspace{-3mu}| |f|^{2r}|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p/(2r)}(T)}^{\delta/r} h^{(1-\delta_0)}
|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p}(T)}^{2(1-\delta)}\\
&
\lesssim h^{1-\delta_0}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p}(T)}^2,
\end{align*}
where $\delta_0=\delta_0(\delta)\to0$ as $\delta\to0$,
and complete the proof.

\end{proof}
Next, we want to give an estimate to $\sup_x\|\mu_s^x-\mu_t^x\|_{var}$. To this end, we will use the relation between the PDE and the SDE. For any $T>0$, consider the following backward PDE:
\begin{align}\label{BPDE0}
\partial_t u^T+B\cdot\nabla u^T+a_{ij}\partial_i\partial_j u^T=0,\quad u^T(T)=\varphi,
\end{align}
where $\varphi\in C^\infty_b$. By Proposition \ref{well-PDE},
there exists a unique solution $u^T$ to \eqref{BPDE0} in the sense of Definition \ref{Pre4:Def}. Denote by
\begin{align*}
P^X_{s,t}f(x):={\mathbb E} f(X^x_{s,t}),\quad P^X_t:=P^X_{0,t}.
\end{align*}
By \eqref{BB06}, the domain of $P^X_{s,t}$ includes $\widetilde{L}^p$ for any $p\in(1,\infty]$.
Then we have the following probability representation;
see \cite[Theorem 6.5.4]{Fried1975} for classical case.
\begin{proposition}\label{Pro:PR}
Let $T>0$, $\varphi\in C^\infty_b({\mathbb R}^d)$, $u^T$ and $X^x_{s,t}$
be the solutions to \eqref{BPDE0} and \eqref{II03} respectively. Then,
\begin{align}\label{PR001}
u^T(s,x)={\mathbb E}\varphi(X^x_{s,T})=P^X_{s,T}\varphi(x).
\end{align}

\end{proposition}
\begin{proof}
It is straightforward to obtain \eqref{PR001} by applying the generalized It\^o's formula \eqref{GIF001}
to function $t\mapsto u^T(t,X^x_{s,t})$ and taking the expectation.
\end{proof}
Apart from the probability representation, by generalized It\^o's formula, we have the following Duhamel's formula.
\begin{lemma}[Duhamel's formula]
For any $\varphi\in C^\infty_b({\mathbb R}^d)$,
\begin{align}\label{DF002}
P^X_{s,t}\varphi(x)=P^\sigma_{t-s}\varphi(x)+\int_s^t P^X_{s,r}\(B(r)\cdot\nabla P^\sigma_{t-r}\varphi\)(x){\mathord{{\rm d}}} r.
\end{align}

\end{lemma}
\begin{proof}
For any $t\in[0,T]$, let $v^t=v^t(t,x)$ be the solution to the following backward PDE:
\begin{align*}
\partial_r v^t+a_{ij}\partial_i\partial_j v^t=0,\quad v^t(t)=\varphi.
\end{align*}
Based on \eqref{PR001}, one sees that $v^t(r)=P^\sigma_{t-r}\varphi$. By generalized It\^o's formula \eqref{GIF001}, we have
\begin{align*}
{\mathbb E} v^t(t,X^x_{s,t})=v^t(s,x)+{\mathbb E}\int_s^t(\partial_rv^t+a_{ij}\partial_i\partial_j v^t+B\cdot\nabla v^t)(r,X^x_{s,r}){\mathord{{\rm d}}} r,
\end{align*}
which implies that
\begin{align*}
P^X_{s,t} \varphi(x)&=P^\sigma_{t-s}\varphi+{\mathbb E}\int_s^t(B(r)\cdot\nabla P^\sigma_{t-r}\varphi)(X^x_{s,r}){\mathord{{\rm d}}} r\\
&=P^\sigma_{t-s}\varphi+\int_s^tP^X_{s,r}(B(r)\cdot\nabla P^\sigma_{t-r}\varphi)(x){\mathord{{\rm d}}} r,
\end{align*}
and we complete the proof.
\end{proof}
It is the position to give the estimate for $\|\mu_t^x-\mu_s^x\|_{var}$.
Notice that similar results have been proved in \cite[Lemma 3.8 (1)(ii)]{Zhao2020}
in the framework of heat kernel estimate. It should be mentioned that
the order of time regularity in \cite{Zhao2020}
depends on the H\"older index $\beta$ of $\sigma$ and is not applicable to our case.
\begin{lemma}\label{NN0000}
Assume \eqref{NBB03}. For any $T>0$ and $q\in[p_0,\infty)$, there is a constant $C=C(\Xi,q)$ such that for any $0<s\leqslant t\leqslant T$ and $\varphi\in C^\infty_b$,
\begin{align}\label{NN000}
\|P^X_t\varphi-P^X_s\varphi\|_\infty\leqslant C\[[(t-s)^{\frac{1}{2}-\frac{d}{2q}} s^{-\frac{1}{2}+\frac{d}{2q}}]\wedge 1\]s^{-\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_q.
\end{align}
In particular, when $q=\infty$, for any $\delta>0$, there is a constant $C=C(\Xi,\delta)$ such that for all $x\in{\mathbb R}^d$,
\begin{align}\label{LL02}
\|\mu^x_t-\mu^x_s\|_{var}\leqslant C\[[(t-s)^{1/2-\delta} s^{-1/2}]\wedge 1\].
\end{align}

\end{lemma}

\begin{proof}
For simplicity, let
\begin{align*}
\alpha:=\frac{1}{2}-\frac{d}{2q}.
\end{align*}
From \eqref{DF002}, one sees that
\begin{align*}
P^X_t\varphi-P^X_s\varphi=&\(P^\sigma_{t}\varphi-P^\sigma_{s}\varphi\)
+\int_s^t P^X_{r}\(B(r)\cdot\nabla P^\sigma_{t-r}\varphi\){\mathord{{\rm d}}} r\\
&+\int_0^s P^X_{r}\Big[B(r)\cdot\nabla \(P^\sigma_{t-r}-P^\sigma_{s-r}\)\varphi\Big]{\mathord{{\rm d}}} r\\
:=&{\mathscr I}_1+{\mathscr I}_2+{\mathscr I}_3.
\end{align*}
Based on \eqref{MM03}, we have
\begin{align*}
\|{\mathscr I}_1\|_\infty &\lesssim \[[(t-s)s^{-1}]\wedge 1\]s^{-\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_q\\
&\lesssim \[(t-s)^\alpha s^{-\alpha}\]s^{-\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_q.
\end{align*}
By \eqref{BB06} and \eqref{BB00}, we have
\begin{align*}
\|{\mathscr I}_2\|_\infty&\lesssim \int_s^t r^{-d/(2p_0)}|\mspace{-3mu}|\mspace{-3mu}| B(r)\cdot\nabla P^\sigma_{t-r}\varphi|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} r\\
&\lesssim \int_s^t r^{-d/(2p_0)}\|\nabla P^\sigma_{t-r}\varphi\|_\infty{\mathord{{\rm d}}} r\\
&\lesssim s^{-d/(2q)}\int_s^t r^{-d/(2p_0)+d/(2q)}(t-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_q\\
&:=s^{-d/(2q)}K(t,s)|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_q,
\end{align*}
where
\begin{align*}
K(t,s)&\leqslant \left(s^{-\frac{d}{2p_0}+\frac{d}{2q}}\int_s^t (t-r)^{\alpha-1}{\mathord{{\rm d}}} r\right)\wedge\left((t-s)^{\frac{1}{2}-\frac{d}{2p_0}}\int_0^1r^{-\frac{d}{2p_0}+\frac{d}{2q}}(1-r)^{-\frac{1}{2}-\frac{d}{2q}}{\mathord{{\rm d}}} r\right)\\
&\lesssim (s^{-\alpha}(t-s)^\alpha)\wedge 1
\end{align*}
provided by $q\geqslant p_0>d$.

The rest one is ${\mathscr I}_3$. By \eqref{BB06} and \eqref{MM03}, we have
\begin{align*}
\|{\mathscr I}_3\|_\infty&\lesssim \int_0^s r^{-d/(2p_0)}|\mspace{-3mu}|\mspace{-3mu}| B(r)\cdot\nabla \(P^\sigma_{t-r}-P^\sigma_{s-r}\)\varphi|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} r\\
&\lesssim \int_0^s r^{-d/(2p_0)}\|\nabla \(P^\sigma_{t-r}-P^\sigma_{s-r}\)\varphi\|_\infty{\mathord{{\rm d}}} r\\
&\lesssim \int_0^s r^{-d/(2p_0)}\[[(t-s)^{\frac12}(s-r)^{-\frac12}]\wedge1\](s-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_q.
\end{align*}
Notice that
\begin{align}
{\mathscr J}&:=\int_0^s r^{-d/(2p_0)}\[[(t-s)^{\frac12}(s-r)^{-\frac12}]\wedge1\](s-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r\nonumber\\
&\leqslant\int_0^{s} r^{-d/(2p_0)}(s-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r\lesssim s^{1/2-d/(2p_0)-d/(2q)}\lesssim s^{-d/(2q)}\label{Star0011}
\end{align}
provided $p_0>d$. In addition, when $r\in(0,\frac{s}{2}]$, one sees that
\begin{align*}
\[(t-s)^{\frac12}(s-r)^{-\frac12}\]\wedge1\leqslant(t-s)^{\frac{1}{2}-\frac{d}{2q}}
(s-r)^{-\frac{1}{2}+\frac{d}{2q}}
\leqslant(t-s)^{\frac{1}{2}-\frac{d}{2q}}\(\frac{s}{2}\)^{-\frac{1}{2}+\frac{d}{2q}}.
\end{align*}
 Hence,
\begin{align*}
{\mathscr J}\lesssim
&
s^{-1/2-d/(2q)}\int_0^{s/2} r^{-d/(2p_0)}\[[(t-s)^{\frac12}(s-r)^{-\frac12}]\wedge1\]{\mathord{{\rm d}}} r\\
&
+ s^{-d/(2p_0)} \int_{s/2}^s \[[(t-s)^{\frac12}(s-r)^{-{\frac12}}]\wedge1\](s-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r\\
\lesssim
&
s^{-1}(t-s)^{1/2-d/(2q)}\int_0^{s/2} r^{-d/(2p_0)}{\mathord{{\rm d}}} r\\
&
+ s^{-d/(2p_0)} \int_{s/2}^s \[[(t-s)^{\frac12}(s-r)^{-\frac12}]\wedge1\](s-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r.
\end{align*}
By a change of variable, we have
 \begin{align*}
{\mathscr J}\lesssim
&
s^{-d/(2p_0)}(t-s)^{1/2-d/(2q)}\\
&
+ s^{-d/(2p_0)}(t-s)^{1-1/2-d/(2q)} \int_{0}^{s/(t-s)} \[(s-r)^{-\frac12}\wedge1\](s-r)^{-1/2-d/(2q)}{\mathord{{\rm d}}} r\\
\lesssim
&
s^{-d/(2p_0)}(t-s)^{1-1/2-d/(2q)}
\lesssim s^{-\frac{1}{2}}(t-s)^{\frac{1}{2}-\frac{d}{2q}}
\end{align*}
provided $q<\infty$ and $p_0>d$, which implies that
\begin{align*}
{\mathscr J}\lesssim \(s^{-\frac12}(t-s)^\alpha\)\wedge s^{-\frac{d}{2q}}=\Big([(t-s)^\alpha s^{-\alpha}]\wedge 1\Big)s^{-\frac{d}{2q}}
\end{align*}
combined with \eqref{Star0011}. Up to now, we obtain \eqref{NN000}.
 In particular, noting that
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| \varphi|\mspace{-3mu}|\mspace{-3mu}|_{q}\lesssim \|\varphi\|_{\infty},\quad \forall q<\infty,
\end{align*}
by Lusin's theorem and \eqref{NN000}, we have for all $\delta>0$ and $q\in[p_0\vee(d/(2\delta)),\infty)$
\begin{align*}
\|\mu^x_t-\mu^x_s\|_{var}&=\sup_{\varphi\in C^\infty_b}\frac{|P^X_t\varphi(x)-P^X_s\varphi(x)|}{\|\varphi\|_\infty}\lesssim [(t-s)^{\frac{1}{2}-\frac{d}{2q}} s^{-\frac{1}{2}+\frac{d}{2q}}]s^{-\frac{d}{2q}}\\
&\lesssim (t-s)^{\frac{1}{2}-\delta} s^{-\frac{1}{2}}.
\end{align*}
Moreover, it is easy to see that
\begin{align*}
\|\mu^x_t-\mu^x_s\|_{var}\leqslant 2,
\end{align*}
which completes the proof.
\end{proof}
\iffalse
 By \cite[Lemma 2.6]{RZ21} and \eqref{Kry02}, we have the following Krylov's estimate in distribution dependent version.
\begin{lemma}
Let $T>0$, $p_1,p_2,q\in(1,\infty)$ with $d/p_1<2$. For any $T>0$ and $Y\in {\mathbb K}^{p_2,q}_{T,\kappa}$, there is a constant $C=C(d,T,|\mspace{-3mu}|\mspace{-3mu}| B|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)},p_1)$ such that for any positive $f\in \widetilde{{\mathbb L}}^{p_1,p_2}_q(T)$,
\begin{align}\label{Kry03}
{\mathbb E}\int_0^T f(r,Y_r,\mu_r){\mathord{{\rm d}}} r\leqslant C\kappa|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_1,p_2}_{q}(T)}.
\end{align}
\end{lemma}
\fi
The following lemma is the distribution dependent version of \eqref{BB05}.
\begin{lemma}
For any $T>0$, $p\in(d\vee 2,\infty)$, assume that $f:{\mathbb R}_+\times{\mathbb R}^{d}\times{\mathcal P}({\mathbb R}^d)\to{\mathbb R}$ satisfying that
\begin{align*}
\kappa_f:=\sup_{t\in[0,T]}\sup_{\mu,\nu}\Big(|\mspace{-3mu}|\mspace{-3mu}| f(t,\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_p+\frac{|\mspace{-3mu}|\mspace{-3mu}| f(t,\cdot,\mu)-f(t,\cdot,\nu)|\mspace{-3mu}|\mspace{-3mu}|_p}{\|\mu-\nu\|_{var}}\Big)<\infty.
\end{align*}
Then, for any $\delta>0$, there is a constant $C=C(\Xi,p,\delta)$ such that for any $x\in{\mathbb R}^d$ and $h\in(0,1)$
\begin{align}\label{EE00}
\begin{split}
{\mathbb E}\Big(\sup_{t\in[0,T]}\Big|\int_{0}^t(f(s,X^x_s,\mu^x_s)&-f(s,X^x_{\pi_h(s)},\mu^x_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\Big)\leqslant C(\kappa_f)^2 h^{1-\delta}.
\end{split}
\end{align}

\end{lemma}
\begin{proof}
For simplicity, we drop the subscript $x$ on $X^x$ and $\mu^x$.
First of all, notice that
\begin{align*}
&~~~{\mathbb E}\Big(\sup_{t\in[0,T]}\Big|\int_{0}^t(f(s,X_s,\mu_s)-f(s,X_{\pi_h(s)},\mu_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\Big)\\
\lesssim&~ {\mathbb E}\Big(\sup_{t\in[0,T]}\Big|\int_{0}^t(f(s,X_s,\mu_s)-f(s,X_{\pi_h(s)},\mu_s)){\mathord{{\rm d}}} s\Big|^2\Big)\\
&+{\mathbb E}\Big(\sup_{t\in[0,T]}\Big|\int_{0}^t(f(s,X_{\pi_h(s)},\mu_s)-f(s,X_{\pi_h(s)},\mu_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\Big)\\
:=&{\mathscr I}^h_1+{\mathscr I}^h_2.
\end{align*}
By \eqref{BB05}, for any $\delta>0$, we have
\begin{align*}
{\mathscr I}^h_1\lesssim h^{1-\delta}\sup_{t\in[0,T]}|\mspace{-3mu}|\mspace{-3mu}| f(t,\cdot,\mu_t)|\mspace{-3mu}|\mspace{-3mu}|_{p}^2\lesssim  h^{1-\delta}(\kappa_f)^2.
\end{align*}
For ${\mathscr I}^h_2$, we use the same method in Step 2 of the proof to Lemma \ref{CruL}. For any $0\leqslant s<t\leqslant T$,
 set $\|\cdot\|:=\({\mathbb E}|\cdot|^2\)^{1/2}$,
\begin{align*}
{\mathscr J}_{s,t}^h:=\int_{s+h}^{t+h}|f(r,X_{\pi_h(r)},\mu_r)-f(r,X_{\pi_h(r)},\mu_{\pi_h(r)})|{\mathord{{\rm d}}} r
\end{align*}
and notice that
\begin{align*}
{\mathscr I}^h_2&\leqslant {\mathbb E}\Big|\int_{h}^{T+h}|f(s,X_{\pi_h(s)},\mu_s)-f(s,X_{\pi_h(s)},\mu_{\pi_h(s)})|{\mathord{{\rm d}}} s\Big|^2=\|{\mathscr J}_{0,T}^h\|^2,
\end{align*}
provided $\pi_h(r)=r$ when $r<h$.
Based on H\"older's inequality and \eqref{Kry02}, one sees that for any $2/q+d/p<1$,
\begin{align*}
{\mathbb E}|{\mathscr J}_{s,t}^h|^2&\lesssim (t-s){\mathbb E}\int_{s+h}^{t+h}
\Big|f(r,X_{\pi_h(r)},\mu_r)-f(r,X_{\pi_h(r)},\mu_{\pi_h(r)})\Big|^2{\mathord{{\rm d}}} r\\
&
\lesssim (t-s)\Big(\int_{s+h}^{t+h}|\mspace{-3mu}|\mspace{-3mu}| f(r,\cdot,\mu_r)-f(r,\cdot,\mu_{\pi_h(r)})|\mspace{-3mu}|\mspace{-3mu}|_{p}^{q}{\mathord{{\rm d}}} r\Big)^{2/q}\\
&
\lesssim (\kappa_f)^2(t-s)\Big(\int_{s+h}^{t+h}\|\mu_r-\mu_{\pi_h(r)}\|_{var}^{q}{\mathord{{\rm d}}} r\Big)^{2/q}.
\end{align*}
Then, by \eqref{LL02} and the fact $q>2$, we have for any $\delta>0$
\begin{align*}
\|{\mathscr J}_{s,t}^h\|&\lesssim \kappa_f(t-s)^{1/2}\Big(\int_{s+h}^{t+h} h^{(\frac12-\delta)q}(\pi_h(r))^{-q/2}{\mathord{{\rm d}}} r\Big)^{1/q}\\
&\lesssim \kappa_f(t-s)^{1/2}h^{1/2-\delta}\Big(\int_{s}^{t} r^{- q/2}{\mathord{{\rm d}}} r\Big)^{1/q}\lesssim \kappa_fh^{1/2-\delta}(t-s)^{1/2+1/q}s^{-1/2}.
\end{align*}
Taking $t_n:=2^{-n}T$, we have
\begin{align*}
\|{\mathscr J}_{0,T}^h\|
&
\leqslant\sum_{n=0}^\infty\|{\mathscr J}_{t_{n+1},t_n}^h\|\lesssim \kappa_fh^{1/2-\delta}
\sum_{n=0}^\infty(t_n-t_{n+1})^{1/2+1/q}{(t_{n+1})}^{-1/2}\\
&
\lesssim \kappa_fh^{1/2-\delta} T^{1/q}\sum_{n=0}^\infty 2^{-\frac{n+1}{q}}\lesssim \kappa_fh^{1/2-\delta},
\end{align*}
which implies ${\mathscr I}^h_2\lesssim (\kappa_f)^2h^{1-2\delta}$
and we complete the proof.

\end{proof}

Now, we have the following strong fluctuation result, which is crucial in this paper.
\begin{lemma}\label{LemEE01}
Let $T>0$ and $g:{\mathbb R}_+\times{\mathbb R}^d\to{\mathbb R}$ be a bounded function satisfying
\begin{align}\label{Cru00}
c_g:=\sup_{t\ne s\in[0,T]\atop x\in{\mathbb R}^d}\Big(|g(t,x)|+\frac{|g(t,x)-g(s,x)|}{|t-s|^\alpha}\Big)<\infty
\end{align}
with some $\alpha>0$.
Assume \eqref{NBB03},  {\bf (H$_b^1$)} and {\bf (H$_b^2$)} hold. For any $\delta>0$, there is a constant $C=C(\Xi,\kappa_0,\delta,c_g)$ such that for any $x\in{\mathbb R}^d$ and $\varepsilon>0$,
\begin{align}\label{EE01}
\begin{split}
&{\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^tg(s,X^x_s)(b(\frac{s}{\varepsilon},X^x_s,\mu^x_s)-\bar{b}(X^x_s,\mu^x_s)){\mathord{{\rm d}}} s\Big|^2\right)\\
&\leqslant C\inf_{h>0}\left(h^{1-\delta}
+h^{2\alpha}+\left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2\right).
\end{split}
\end{align}

\end{lemma}

\begin{proof}
For simplicity, we drop the subscript $x$ on $X^x$ and $\mu^x$.
Set $b_\varepsilon(t):=b(t/\varepsilon)$ and
\begin{align*}
{\bf X}_\cdot:=(X_\cdot,\mu_\cdot).
\end{align*}
For any $h\in(0,1)$,
\begin{align*}
&{\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^tg(s,X^x_s)(b_\varepsilon(s,X^x_s,\mu^x_s)-\bar{b}(X^x_s,\mu^x_s)){\mathord{{\rm d}}} s\Big|^2\right)\\
&
\lesssim {\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^t((gb_\varepsilon)(s,{\bf X}_s)-(gb_\varepsilon)(s,{\bf X}_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)\\
&\quad
+{\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^tg(s,X_{\pi_h(s)})(b_\varepsilon(s,{\bf X}_{\pi_h(s)})-\bar{b}({\bf X}_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)\\
&\quad
+{\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^t((g\bar{b})(s,{\bf X}_s)-(g\bar{b})(s,{\bf X}_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)\\
&
:={\mathscr I}^{\varepsilon,h}_1+{\mathscr I}^{\varepsilon,h}_2+{\mathscr I}^{\varepsilon,h}_3.
\end{align*}
Based on condition {\bf (H$^1_b$)}, by \eqref{EE00}, for any $\delta>0$, we have
\begin{align*}
{\mathscr I}^{\varepsilon,h}_1+{\mathscr I}^{\varepsilon,h}_3
\lesssim \|g\|_\infty^2\kappa_0^2h^{1-\delta}.
\end{align*}

For ${\mathscr I}^{\varepsilon,h}_2$, notice that by \eqref{Cru00} and \eqref{Kry02} with $p_0>(d/2)\vee1$,
\begin{align*}
{\mathscr I}^{\varepsilon,h}_2&\lesssim {\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^tg(\pi_h(s),X_{\pi_h(s)})(b_\varepsilon(s,{\bf X}_{\pi_h(s)}\)-\bar{b}({\bf X}_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)\\
&+{\mathbb E}\left(\sup_{t\in[0,T]}\int_0^t|s-\pi_h(s)|^{2\alpha}\Big(|b_\varepsilon(s,{\bf X}_{\pi_h(s)})|^2+|\bar{b}({\bf X}_{\pi_h(s)}))|^2\Big){\mathord{{\rm d}}} s\right)\\
&\lesssim {\mathscr I}^{\varepsilon,h}_{21}+h^{2\alpha},
\end{align*}
where
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{21}:={\mathbb E}\left(\sup_{t\in[0,T]}\Big|\int_0^tg(\pi_h(s),X_{\pi_h(s)})(b_\varepsilon(s,{\bf X}_{\pi_h(s)})-\bar{b}({\bf X}_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right).
\end{align*}
It suffices to show
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{21}\lesssim h+\(\omega\left(\frac{h}{\varepsilon}\right)\)^2
\end{align*}
in what follows.
Letting $M(t)=[t/h]$ and noting that $\pi_h(s)=s$ for $s\in[0,h)$, we have
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{21}
\lesssim
&
{\mathbb E}\left(\sup_{t\in[0,h)}\Big|\int_0^tg(s,X_s)(b_\varepsilon(s,{\bf X}_{s})-\bar{b}({\bf X}_{s})){\mathord{{\rm d}}} s\Big|^2\right)\\
&
+{\mathbb E}\left(\sup_{t\in[h,T]}\Big|\int_0^hg(s,X_s)(b_\varepsilon(s,{\bf X}_{s})-\bar{b}({\bf X}_{s})){\mathord{{\rm d}}} s\Big|^2\right)\\
&
+{\mathbb E}\left(\sup_{t\in[h,T]}\Big|\int_{M(t)h}^tg(\pi_h(s),X_{\pi_h(s)})(b_\varepsilon(s,{\bf X}_{\pi_h(s)})
-\bar{b}({\bf X}_{\pi_h(s)}){\mathord{{\rm d}}} s\Big|^2\right)\\
&
+{\mathbb E}\left(\sup_{t\in[h,T]}\Big|\int_h^{M(t)h}g(\pi_h(s),X_{\pi_h(s)})(b_\varepsilon(s,{\bf X}_{\pi_h(s)})
-\bar{b}({\bf X}_{\pi_h(s)})){\mathord{{\rm d}}} s\Big|^2\right)\\
:=&
{\mathscr I}^{\varepsilon,h}_{211}+{\mathscr I}^{\varepsilon,h}_{212}+{\mathscr I}^{\varepsilon,h}_{213}+{\mathscr I}^{\varepsilon,h}_{214}.
\end{align*}
It follows from H\"older's inequality
 that
\begin{align*}
\sum_{i=1}^3{\mathscr I}^{\varepsilon,h}_{21i}
&
\lesssim h \|g\|_\infty^2\left({\mathbb E}\left[\int_0^T|b_\varepsilon(s,{\bf X}_{s})|^2+|\bar{b}({\bf X}_{s})|^2{\mathord{{\rm d}}} s\right]
+{\mathbb E}\left[\int_0^T|b_\varepsilon(s,{\bf X}_{\pi_h(s)})|^2+|\bar{b}({\bf X}_{\pi_h(s)})|^2{\mathord{{\rm d}}} s\right]\right)\\
&
\lesssim h\|g\|_\infty^2\kappa_0^2,
\end{align*}
provided by \eqref{Kry02}.
Thus, we only need to prove
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{124}\lesssim \(\omega\left(\frac{h}{\varepsilon}\right)\)^2.
\end{align*}
By the definition of $\pi_h$, it is easy to see that
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{214}\leqslant{\mathbb E}\Big(\sup_{2\leqslant m\leqslant M(T)}\Big|\sum_{k=1}^{m-1}g(kh,X_{kh})\int_{kh}^{(k+1)h}(b_\varepsilon(s,{\bf X}_{kh})-\bar{b}({\bf X}_{kh})){\mathord{{\rm d}}} s\Big|^2\Big).
\end{align*}

Based on the fact $|\sum_{k=1}^{m-1}a_k|^2\leqslant (m-1)\sum_{k=1}^{m-1}|a_k|^2$, one sees that
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{214}
&
\leqslant M(T)c_g^2\sum_{k=1}^{M(T)-1}{\mathbb E}\Big|\int_{kh}^{(k+1)h}
(b_\varepsilon(s,{\bf X}_{kh})-\bar{b}({\bf X}_{kh})){\mathord{{\rm d}}} s\Big|^2.
\end{align*}
By a change of variable and \eqref{in:00}, we have
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{214}
&
\lesssim
M(T)\sum_{k=1}^{M(T)-1}{\mathbb E}\left|\varepsilon\int_{kh/\varepsilon}^{(k+1)h/\varepsilon}\(b(s,{\bf X}_{kh})-\bar{b}({\bf X}_{kh})\){\mathord{{\rm d}}} s\right|^2\\
&
\lesssim [\frac{T}{h}]h^2\sum_{k=1}^{M(T)-1}{\mathbb E}\Big|\frac{\varepsilon}{h}\int_{kh/\varepsilon}^{(k+1)h/\varepsilon}
(b(s,{\bf X}_{kh})-\bar{b}({\bf X}_{kh})){\mathord{{\rm d}}} s\Big|^2\\
&
\lesssim h \left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2
\sum_{k=1}^{M(T)-1}{\mathbb E} |H({\bf X}_{kh})|^2.
\end{align*}
Notice that
\begin{align*}
h\sum_{k=1}^{M(T)-1}{\mathbb E} |H({\bf X}_{kh})|^2={\mathbb E}\int_{h}^{M(T)h}|H({\bf X}_{\pi_h(s)})|^2{\mathord{{\rm d}}} s
\leqslant{\mathbb E}\int_{h}^{T}|H({\bf X}_{\pi_h(s)})|^2{\mathord{{\rm d}}} s.
\end{align*}
By \eqref{Kry02} again, we have
\begin{align*}
{\mathscr I}^{\varepsilon,h}_{214}\lesssim \left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2
\sup_{\mu}|\mspace{-3mu}|\mspace{-3mu}| H(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}^2
\end{align*}
and complete the proof.

\end{proof}

\subsection{Time regularity for solutions to parabolic equation}\label{SubS32}
In this section, we establish the time regularity for the solution of PDE \eqref{NBB00}.
Firstly, we give the following probability representation of the solution to PDE when $B\equiv0$.
\begin{lemma}
Let $B\equiv0$ and $u$ be a solution to PDE \eqref{NBB00}. Then,
\begin{align}\label{CCA01}
u(t)=\int_0^t e^{-\lambda (t-s)}P^\sigma_{t-s}f(s){\mathord{{\rm d}}} s
+e^{-\lambda t}P^\sigma_t\varphi.
\end{align}
\end{lemma}
\begin{proof}
Using the generalized It\^o's formula \eqref{GIF001} to
$$
s\to e^{-\lambda s}u(t-s,Z_s^x),
$$
we get
\begin{align*}
e^{-\lambda t}u(0,Z_t^x)-u(t,x)=&\int_0^t e^{-\lambda s}(-\partial_su+a_{ij}\partial_i\partial_j u-\lambda u)(t-s,Z_s^x){\mathord{{\rm d}}} s\\
&+\int_0^t e^{-\lambda s}\nabla u(t-s,Z_s^x){\mathord{{\rm d}}} W_s.
\end{align*}
Taking the expectations of both sides, we obtain that
\begin{align*}
e^{-\lambda t}P^\sigma_t\varphi-u(t)=-\int_0^t e^{-\lambda s}P^\sigma_sf(t-s){\mathord{{\rm d}}} s,
\end{align*}
which is \eqref{CCA01} by a change of variable and this completes the proof.
\end{proof}

Using the above Lemma, we have the following time H\"older regularity of $\nabla u$.
\begin{lemma}
Assume $\varphi\equiv0$. Under the condition \eqref{NBB03} with some $p_0\in(d,\infty)$, for any $\lambda\ge0$, there is a constant $C=C(\Xi,p,\lambda)$ such that for all $t,s\in[0,T]$, $f\in\widetilde{{\mathbb L}}^{p_0}(T)$ the solution $u$ to \eqref{NBB00} in the sense of Definition \ref{Pre4:Def} satisfying
\begin{align}\label{HH02}
\|\nabla u(t)-\nabla u(s)\|_\infty\leqslant C|t-s|^{\frac{1}{2}-\frac{d}{2p_0}}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}.
\end{align}
\end{lemma}
\begin{remark}\label{Re002}\rm
By Remark \ref{Re001} and \eqref{HH02}, we further have there is a version of solution such that $u\in C([0,T];{\mathbf C}^1)$.
\end{remark}
\begin{proof}
Firstly, since $B,f\in\widetilde{{\mathbb L}}^{p_0}(T)\subset\widetilde{{\mathbb L}}^{p_0}_q(T)$, $\forall q$, we indeed have the unique solution $u$.
Set
\begin{align*}
g(s):=B\cdot\nabla u(s)+f(s).
\end{align*}
In view of \eqref{CCA022},
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}| g|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}
\leqslant|\mspace{-3mu}|\mspace{-3mu}| b|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}\|\nabla u\|_{{\mathbb L}^\infty_T}
+|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}\lesssim |\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}.
\end{align*}
Then, for any $0\leqslant s<t\leqslant T$,
by \eqref{CCA01}, one sees that
\begin{align*}
\|\nabla u(t)-\nabla u(s)\|_\infty
&
=\left\|\int_s^t e^{-\lambda(t-r)}\nabla P^\sigma_{t-r}g(r){\mathord{{\rm d}}} r\right\|_\infty\\
&
+\left\|\int_0^s \(e^{-\lambda(t-r)}-e^{-\lambda(s-r)}\)\nabla P^\sigma_{t-r}g(r){\mathord{{\rm d}}} r\right\|_\infty\\
&
+\left\|\int_0^s e^{-\lambda(s-r)}\(\nabla P^\sigma_{t-r}-\nabla P^\sigma_{s-r}\)g(r){\mathord{{\rm d}}} r\right\|_\infty\\
&
=:{\mathscr I}_1+{\mathscr I}_2+{\mathscr I}_3.
\end{align*}
By \eqref{BB00}, one sees that
\begin{align*}
{\mathscr I}_1
&
\lesssim \int_s^t (t-r)^{-1/2-d/(2{p_0})}|\mspace{-3mu}|\mspace{-3mu}| g(r)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} r\nonumber\\
&
\lesssim \int_s^t (t-r)^{-1/2-d/(2{p_0})}{\mathord{{\rm d}}} r|\mspace{-3mu}|\mspace{-3mu}| g|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}\nonumber\\
&
\lesssim|t-s|^{\frac{1}{2}-\frac{d}{2{p_0}}}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}
\end{align*}
For ${\mathscr I}_2$, noting that $|e^{-x}-e^{-y}|\leqslant |x-y|$ for any $x,y>0$, it follows from \eqref{BB00} that
\begin{align*}
{\mathscr I}_2
&
\lesssim |t-s|\int_0^s (t-r)^{-1/2-d/(2{p_0})}|\mspace{-3mu}|\mspace{-3mu}| g(r)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} r\\
&
\lesssim|t-s||\mspace{-3mu}|\mspace{-3mu}| g|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}\lesssim|t-s||\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}.
\end{align*}
For ${\mathscr I}_3$, by \eqref{MM03} with $q=\infty$, we have
\begin{align*}
{\mathscr I}_3
&
\lesssim \int_0^s([(t-s)^{\frac12}(s-r)^{-\frac12}]\wedge1)(s-r)^{-1/2-d/(2p_0)}{\mathord{{\rm d}}} r|\mspace{-3mu}|\mspace{-3mu}| g|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}.
\end{align*}
By a change of variable, we have
\begin{align*}
{\mathscr I}_3&\lesssim |t-s|^{\frac12-\frac{d}{2p_0}} \int_0^\infty ([(s-r)^{-\frac12}]\wedge1)(s-r)^{-1/2-d/(2p_0)}{\mathord{{\rm d}}} r|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}\\
&\lesssim |t-s|^{1/2-d/(2p_0)}|\mspace{-3mu}|\mspace{-3mu}| f|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}
\end{align*}
and complete the proof.

\end{proof}
Moreover, we also have an estimate of time regularity for solution to the following Cauchy problem.
\begin{lemma}\label{Lemma3.12}
Assume \eqref{NBB03}. Let $\varphi\in C^\infty_b$ and $u$ be the unique solution to the following Cauchy equation on $[0,T]$ in the sense of Definition \ref{Pre4:Def}
\begin{align}\label{CE00}
\partial_t u=a_{ij}\partial_i\partial_j u+B\cdot\nabla u,\quad u_0=\varphi.
\end{align}
Then there is a constant $C=C(\Xi)$ such that for all $0\leqslant s<t\leqslant T$,
\begin{align}\label{NN01}
\|\nabla u(t)-\nabla u(s)\|_\infty\leqslant C (t-s)^{\frac12-\frac{d}{2p_0}}s^{-1+\frac{d}{2p_0}}\|\varphi\|_\infty.
\end{align}

\end{lemma}
\begin{proof}
Firstly, by \eqref{CCA01}, we have
\begin{align}\label{MM07}
u(t)=\int_0^t P^\sigma_{t-s}(B\cdot \nabla u)(s){\mathord{{\rm d}}} s+P^\sigma_t \varphi,
\end{align}
which implies that
\begin{align*}
\|\nabla u(t)\|_\infty
&
\lesssim \int_0^t (t-s)^{-\frac{1}{2}-\frac{d}{2p_0}}|\mspace{-3mu}|\mspace{-3mu}| B\cdot \nabla u(s)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} s+t^{-\frac{1}{2}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_\infty\\
&
\lesssim \int_0^t (t-s)^{-\frac{1}{2}-\frac{d}{2p_0}}\|\nabla u(s)\|_\infty{\mathord{{\rm d}}} s+t^{-\frac{1}{2}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_\infty
\end{align*}
provided by \eqref{BB00}.
Hence, by Gronwall's inequality, we have
\begin{align}\label{MM06}
\|\nabla u(t)\|_\infty\lesssim t^{-\frac{1}{2}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_\infty.
\end{align}
 \iffalse
From \eqref{MM07}, for any $0\leqslant s<t\leqslant T$,
\begin{align*}
\|u(t)-u(s)\|_\infty\leqslant& \int_s^t\|P_{t-r} (b\cdot\nabla u)(r)\|_\infty{\mathord{{\rm d}}} r+\int_0^s\|(P_{t-r}-P_{s-r})b\cdot\nabla u(r)\|_\infty{\mathord{{\rm d}}} r\\
&+\|(P_t-P_s)\varphi\|_\infty.
\end{align*}
It follows from \eqref{BB00}, \eqref{MM03} and \eqref{MM06} that
\begin{align*}
\|u(t)-u(s)\|_\infty\lesssim& \int_s^t(t-r)^{-\frac{d}{2p}}\|\nabla u(r)\|_\infty{\mathord{{\rm d}}} r\\
&
+\int_0^s\Big\{\[(t-s)(s-r)^{-1}\]\wedge1\Big\}(s-r)^{-\frac{d}{2p}}\|\nabla u(r)\|_\infty{\mathord{{\rm d}}} r\\
&
+\Big\{\[(t-s)s^{-1}\]\wedge1\Big\}s^{-\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_q\\
\lesssim
&
\Big(\int_s^t(t-r)^{-\frac{d}{2p}}r^{-\frac{1}{2}-\frac{d}{2q}}{\mathord{{\rm d}}} r+(t-s)^{\frac{1}{2}}\int_0^s(s-r)^{-\frac{1}{2}-\frac{d}{2p}}r^{-\frac{1}{2}-\frac{d}{2q}}{\mathord{{\rm d}}} r\\
&
+(t-s)^{\frac{1}{2}}s^{-\frac{1}{2}-\frac{d}{2q}}\Big)|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_q\\
\lesssim
& (t-s)^{\frac{1}{2}}s^{-\frac{1}{2}-\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_q,
\end{align*}
which is \eqref{NN00} with the following observation
\begin{align*}
\|u(t)-u(s)\|_\infty\leqslant\|u(t)\|_\infty+\|u(s)\|_\infty\lesssim(t^{-\frac{d}{2q}}+s^{-\frac{d}{2q}})|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_q\lesssim s^{-\frac{d}{2q}}|\mspace{-3mu}|\mspace{-3mu}|\varphi|\mspace{-3mu}|\mspace{-3mu}|_q.
\end{align*}
\fi
By \eqref{MM07}, \eqref{BB00} and \eqref{MM03}, one sees that
\begin{align*}
\|\nabla u(t)-\nabla u(s)\|_\infty\leqslant& \int_s^t\|\nabla P^\sigma_{t-r} (B\cdot\nabla u)(r)\|_\infty{\mathord{{\rm d}}} r+\int_0^s\|\nabla(P^\sigma_{t-r}-P^\sigma_{s-r})B\cdot\nabla u(r)\|_\infty{\mathord{{\rm d}}} r\\
&+\|\nabla(P^\sigma_t-P^\sigma_s)\varphi\|_\infty.\\
\lesssim& \int_s^t(t-r)^{-\frac12-\frac{d}{2p_0}}\|\nabla u(r)\|_\infty{\mathord{{\rm d}}} r\\
&+\int_0^s\Big\{\[(t-s)^{\frac12}(s-r)^{-\frac12}\]\wedge1\Big\}(s-r)^{-\frac12-\frac{d}{2p_0}}\|\nabla u(r)\|_\infty{\mathord{{\rm d}}} r\\
&+\Big\{\[(t-s)^{\frac12}s^{-\frac12}\]\wedge1\Big\}s^{-\frac{1}{2}}\|\varphi\|_\infty\\
:=&{\mathscr I}_1+{\mathscr I}_2+{\mathscr I}_3.
\end{align*}
Then, noting that $x\wedge1\leqslant x^\theta$ for all $x\ge0$ and $\theta\in[0,1]$, we have
\begin{align*}
{\mathscr I}_3\lesssim \left\{\[(t-s)^{\frac12}s^{-\frac12}\]\wedge1\right\}s^{-\frac{1}{2}}\|\varphi\|_\infty
\lesssim (t-s)^{\frac12-\frac{d}{2p_0}}s^{-1+\frac{d}{2p_0}}\|\varphi\|_\infty.
\end{align*}

Based on \eqref{MM06} and a change of variable,
\begin{align*}
{\mathscr I}_1&\lesssim \|\varphi\|_\infty\int_s^t(t-r)^{-\frac12-\frac{d}{2p_0}}r^{-\frac12}{\mathord{{\rm d}}} r\lesssim (t-s)^{\frac12-\frac{d}{2p_0}}s^{-\frac12}\|\varphi\|_\infty\\
&\lesssim (t-s)^{\frac12-\frac{d}{2p_0}}s^{-1+\frac{d}{2p_0}}\|\varphi\|_\infty
\end{align*}
provided $p_0>d$.
For ${\mathscr I}_2$, by \eqref{MM06} again, we divide $[0,s]$ into $[0,s/2)$ and $[s/2,s]$ and have
\begin{align*}
\frac{{\mathscr I}_2}{\|\varphi\|_\infty}&\lesssim\(\int_0^{s/2}+\int_{s/2}^s\)\Big\{\[(t-s)^{\frac12}(s-r)^{-\frac12}\]\wedge1\Big\}(s-r)^{-\frac12-\frac{d}{2p_0}}r^{-\frac12}{\mathord{{\rm d}}} r\\
&\lesssim (t-s)^{\frac12-\frac{d}{2p_0}}s^{-1}\int_0^{s/2} r^{-\frac12}{\mathord{{\rm d}}} r+s^{-\frac12}\int_{s/2}^s \Big\{\[(t-s)^{\frac12}(s-r)^{-\frac12}\]\wedge1\Big\}(s-r)^{-\frac12-\frac{d}{2p_0}}{\mathord{{\rm d}}} r\\
&\lesssim (t-s)^{\frac12-\frac{d}{2p_0}}s^{-\frac12}+s^{-\frac12}\int_{0}^s \Big\{\[(t-s)^{\frac12}r^{-\frac12}\]\wedge1\Big\}r^{-\frac12-\frac{d}{2p_0}}{\mathord{{\rm d}}} r,
\end{align*}
where we used the fact
\begin{align*}
&\left[\Big((t-s)^{\frac12}(s-r)^{-\frac12}\Big)\wedge1\right](s-r)^{-\frac12-\frac{d}{2p_0}}\\
&\leqslant (t-s)^{\frac12-\frac{d}{2p_0}}(s-r)^{-1}\leqslant 2 (t-s)^{\frac12-\frac{d}{2p_0}}s^{-1}\quad \forall r\in[0,s/2).
\end{align*}
 From a change of variable, we have
\begin{align*}
{\mathscr I}_2
&\lesssim \Big[(t-s)^{\frac12-\frac{d}{2p_0}}s^{-\frac12}
+s^{-\frac12}(t-s)^{1-\frac12-\frac{d}{2p_0}}\int_{0}^\infty (r^{-\frac12}\wedge1)r^{-\frac12-\frac{d}{2p_0}}{\mathord{{\rm d}}} r\Big]\|\varphi\|_\infty\\
&\lesssim (t-s)^{\frac12-\frac{d}{2p_0}}s^{-\frac12}\|\varphi\|_\infty\lesssim(t-s)^{\frac12-\frac{d}{2p_0}}s^{-1+\frac{d}{2p_0}}\|\varphi\|_\infty
\end{align*}
provided $p_0\in(d,\infty)$ and complete the proof.

\end{proof}

\section{Convergence rate of total variation}\label{Section4}
In this section, under {\bf (H$_\sigma$)}, {\bf (H$_b^1$)} and {\bf (H$_b^2$)}, we will give the convergence rate of $\|\mu^\varepsilon-\mu\|_{var}$.
Recall that in probability space $(\Omega,{\mathscr F},{\mathbb P},({\mathscr F}_s)_{s\ge0})$ we have the unique strong solution $(X^\varepsilon_\cdot,X_\cdot)$ to the following systems
\begin{align}\label{ori:DDSDE}
{\mathord{{\rm d}}} X^\varepsilon_t=b(t/\varepsilon,X^\varepsilon_t,\mu^\varepsilon_t){\mathord{{\rm d}}} t+\sigma(X^\varepsilon_t){\mathord{{\rm d}}} W_t, \ \ X^\varepsilon_0=\xi,
\end{align}
and
\begin{align}\label{ave:DDSDE}
{\mathord{{\rm d}}} X_t=\bar{b}(X_t,\mu_t){\mathord{{\rm d}}} t+\sigma(X_t){\mathord{{\rm d}}} W_t, \ \ X_0=\xi,
\end{align}
where $\mu^\varepsilon_t$ and $\mu_t$ are the distributions of $X^\varepsilon_t$ and $X_t$ respectively.
For simplicity, in the sequel, let
\begin{align*}
 b_\varepsilon(t):=b(t/\varepsilon)~~(\varepsilon>0), \ \ \text{and}\ \ b_0:=\bar{b}.
\end{align*}
For any $x\in{\mathbb R}^d$ and $t\geqslant s\geq0$,
let $(Y^\varepsilon_{s,t}(x),Y_{s,t}(x))$ be the unique
strong solution to the following SDEs
\begin{align*}
{\mathord{{\rm d}}} Y^\varepsilon_{s,t}(x)=b_\varepsilon(t,Y^\varepsilon_{s,t}(x),\mu^\varepsilon_t){\mathord{{\rm d}}} t+\sigma(Y^\varepsilon_{s,t}(x)){\mathord{{\rm d}}} W_t,\quad Y^\varepsilon_{s,s}(x)=x
\end{align*}
and
\begin{align*}
{\mathord{{\rm d}}} Y_{s,t}(x)=\bar{b}(Y_{s,t}(x),\mu_t){\mathord{{\rm d}}} t+\sigma(Y_{s,t}(x)){\mathord{{\rm d}}} W_t,\quad Y_{s,s}(x)=x.
\end{align*}
Set $Y^\varepsilon_{t}(x):=Y^\varepsilon_{0,t}(x)$ and $Y_{t}(x):=Y_{0,t}(x)$ for all $t\geq0$ and $x\in\R^d$. Denote by $P^{x,\varepsilon}$ and $P^x$ the distribution of $Y^{\varepsilon}_\cdot(x)$ and $Y_\cdot(x)$ in $C([0,T];{\mathbb R}^d)$ respectively. Based on the strong uniqueness of above SDEs, we have
\begin{align}\label{DStar00}
\int_{{\mathbb R}^d} P^{x,\varepsilon}~~{\mathbb P}\circ\xi^{-1}({\mathord{{\rm d}}} x)={\mathbb P}\circ (X^\varepsilon_\cdot)^{-1}\quad \text{and}\quad\int_{{\mathbb R}^d} P^{x}~~{\mathbb P}\circ\xi^{-1}({\mathord{{\rm d}}} x)={\mathbb P}\circ (X_\cdot)^{-1}.
\end{align}
 Therefore, the estimates in Section \ref{Sec3.1} hold for $X^\varepsilon$ and $X$ where the constants are independent of $\varepsilon$ since
\begin{align*}
\sup_{\varepsilon\ge0}\sup_{\mu\in{\mathcal P}({\mathbb R}^d)}|\mspace{-3mu}|\mspace{-3mu}| b_\varepsilon(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}_T}<\infty.
\end{align*}
Moreover, for any $t\in{\mathbb R}_+$ and $\varphi\in C^\infty_b$, consider the following Kolmogorov's backward equation
\begin{align}\label{BKeq}
\partial_s u^t+a_{ij}\partial_i\partial_j u^t+b_0(\cdot,\mu_s)\cdot\nabla u^t=0,
\end{align}
with the end condition
$$u^t(t)=\varphi.$$
By Proposition \ref{well-PDE} and \ref{Pro:PR}, there exists a unique solution $u^t$ to
\eqref{BKeq} which is given by
$$
u^t(s,x)={\mathbb E}\varphi(Y_{s,t}(x)).
$$
Define $\tilde{u}(s,x)=u^t(t-s,x)$. Then $\tilde{u}$ is the solution to
\begin{equation*}
  \partial_s\tilde{u}=a_{ij}\partial_i\partial_j \tilde{u}
  +b_0(\cdot,\mu_{t-s})\cdot\nabla \tilde{u}, \quad \tilde{u}_0=\varphi.
\end{equation*}
By Lemma \ref{Lemma3.12}, we have
\begin{align}\label{NN04}
\|\nabla u^t(s)\|_\infty \lesssim (t-s)^{-\frac12}\|\varphi\|_\infty
\end{align}
and
\begin{align}\label{FStar00}
\|\nabla u^t(s_1)-\nabla u^t(s_2)\|_\infty\lesssim |s_1-s_2|^{\frac12-\frac{d}{2p_0}}(s_1\wedge s_2)^{-1+\frac{d}{2p_0}}\|\varphi\|_\infty.
\end{align}
Then, for any $t\in[0,T]$, by applying the generalized It\^o's formula \eqref{GIF001} to $u^t(s,Y^\varepsilon_s(x))$, one sees that
\begin{align}\label{NCZ00}
{\mathbb E} \varphi(Y^\varepsilon_t(x))-u^t(0,x)
={\mathbb E}\int_0^t\Big(b_\varepsilon(s,Y^\varepsilon_s(x),\mu^\varepsilon_s)
-b_0(Y^\varepsilon_s(x),\mu_s)\Big)\cdot\nabla u^t(s,Y^\varepsilon_s(x)){\mathord{{\rm d}}} s.
\end{align}
Noting that $u^t(0,x)={\mathbb E} \varphi(Y_t(x))$, based on \eqref{DStar00}, we have
\begin{align}\label{NN03}
\begin{split}
|{\mathbb E} \varphi(X^\varepsilon_t)-{\mathbb E} \varphi(X_t)|&=\Big|\int_{{\mathbb R}^d}\({\mathbb E} \varphi(Y^\varepsilon_t(x))-{\mathbb E} \varphi(Y_t(x))\){\mathbb P}\circ\xi^{-1}({\mathord{{\rm d}}} x)\Big|\\
&\leqslant\sup_x\Big|{\mathbb E}\int_0^t\Big(b_\varepsilon(s,Y^\varepsilon_s(x),\mu^\varepsilon_s)-b_0(Y^\varepsilon_s(x),\mu_s)\Big)\cdot\nabla u^t(s,Y^\varepsilon_s(x)){\mathord{{\rm d}}} s\Big|.
\end{split}
\end{align}

Here is the main result in this section:
\begin{theorem}\label{S4:Main}
Under the conditions {\bf (H$_\sigma$)} and {\bf (H$_b^1$)}--{\bf (H$_b^2$)}. For any $T>0$,
there is a constant $C=C(\kappa_0,\kappa_1,d,T,p_0,\beta)>0$
such that for all $\varepsilon>0$ and $t\in[0,T]$,
\begin{align}\label{LL00}
\|\mu^\varepsilon_t-\mu_t\|_{var}\leqslant C \inf_{h>0}\Big(h^{\frac12-\frac{d}{2p_0}}+\omega\left(\frac{h}{\varepsilon}\right)\Big).
\end{align}
\end{theorem}
\begin{proof}
For the simplicity, in the whole proof, we assume $\|\varphi\|_\infty=1$ and drop the subscript $t$ on $u^t$.
Firstly, let
\begin{align*}
{\mathcal B}^\varepsilon:=\Big|{\mathbb E}\int_0^t\Big(b_\varepsilon(s,Y^\varepsilon_s(x),\mu^\varepsilon_s)-b_0(Y^\varepsilon_s(x),\mu_s)\Big)\cdot\nabla u(s,Y^\varepsilon_s(x)){\mathord{{\rm d}}} s\Big|
\end{align*}
and
\begin{align*}
{\mathcal E}^\varepsilon_h:=\Big|{\mathbb E}\int_0^t\Big(b_\varepsilon(s,Y^\varepsilon_{\pi_h(s)}(x),\mu^\varepsilon_{\pi_h(s)})-b_0(Y^\varepsilon_{\pi_h(s)}(x),\mu^\varepsilon_{\pi_h(s)})\Big)\cdot\nabla u({\pi_h(s)},Y^\varepsilon_{\pi_h(s)}(x)){\mathord{{\rm d}}} s\Big|,
\end{align*}
where $h\in(0,1)$.
For any field $f:{\mathbb R}_+\times{\mathbb R}^d\times{\mathcal P}({\mathbb R}^d)\to{\mathbb R}^d$ and $h\in(0,1)$, define
\begin{align*}
&U_{1,h}^\varepsilon(f):=\Big|{\mathbb E}\int_0^t\((f\cdot\nabla u)(s,Y^\varepsilon_s(x),\mu^\varepsilon_s)-(f\cdot\nabla u)(s,Y^\varepsilon_{\pi_h(s)}(x),\mu^\varepsilon_s)\){\mathord{{\rm d}}} s\Big|\\
&U_{2,h}^\varepsilon(f):=\Big|{\mathbb E}\int_0^t\Big(
\[f(s,Y^\varepsilon_{\pi_h(s)}(x),\mu^\varepsilon_s)
-f(s,Y^\varepsilon_{\pi_h(s)}(x),\mu^\varepsilon_{\pi_h(s)})\]\cdot\nabla u(s,Y^\varepsilon_{\pi_h(s)}(x))\Big){\mathord{{\rm d}}} s\Big|\\
&U_{3,h}^\varepsilon(f):=\Big|{\mathbb E}\int_0^t\Big(
f(s,Y^\varepsilon_{\pi_h(s)}(x),\mu^\varepsilon_{\pi_h(s)})\cdot\[\nabla u(s,Y^\varepsilon_{\pi_h(s)}(x))-\nabla u(\pi_h(s),Y^\varepsilon_{\pi_h(s)}(x))\]\Big){\mathord{{\rm d}}} s\Big|.
\end{align*}
Then, we have
\begin{align*}
{\mathcal B}^\varepsilon\leqslant&{\mathcal E}^\varepsilon_h+\sum_{i=1}^3 \[U_{i,h}^\varepsilon(b_\varepsilon)+U_{i,h}^\varepsilon(b_0)\]\\
&
+\Big|{\mathbb E}\int_0^t\Big(b_0(Y^\varepsilon_s(x),\mu^\varepsilon_s)
-b_0(Y^\varepsilon_s(x),\mu_s)\Big)\cdot\nabla u(s,Y^\varepsilon_s(x)){\mathord{{\rm d}}} s\Big|.
\end{align*}
It follows from \eqref{NN04}, \eqref{BB06} and {\bf (H$_b^1$)} that
\begin{align*}
&\quad\Big|{\mathbb E}\int_0^t\Big(b_0(Y^\varepsilon_s(x),\mu^\varepsilon_s)
-b_0(Y^\varepsilon_s(x),\mu_s)\Big)\cdot\nabla u(s,Y^\varepsilon_s(x)){\mathord{{\rm d}}} s\Big|\\
&\lesssim \int_0^ts^{-\frac{d}{2p_0}}(t-s)^{-\frac12}\|\mu^\varepsilon_s-\mu_s\|_{var}{\mathord{{\rm d}}} s,
\end{align*}
which implies that
\begin{align}\label{NN05}
{\mathcal B}^\varepsilon\lesssim  {\mathcal E}^\varepsilon_h+\sum_{i=1}^3 \[U_{i,h}^\varepsilon(b_\varepsilon)+U_{i,h}^\varepsilon(b_0)\]+\int_0^ts^{-\frac{d}{2p_0}}(t-s)^{-\frac12}\|\mu^\varepsilon_s-\mu_s\|_{var}{\mathord{{\rm d}}} s.
\end{align}

Now, we divide the rest of the proof into two steps. In step 1, we estimate $U_{i,h}^\varepsilon(b_\varepsilon)+U_{i,h}^\varepsilon(b_0)$, $i=1,2,3$, one by one; In step 2, we calculate ${\mathcal E}^\varepsilon_h$ under the assumption {\bf (H$_b^2$)}.

{\bf(Step 1)} We only estimate $U_{i,h}^\varepsilon(b_\varepsilon)$, $U_{i,h}^\varepsilon(b_0)$ is totally the same. Firstly, we estimate $U_{1,h}^\varepsilon(b_\varepsilon)$. By \eqref{NN000} and \eqref{NN04}, we have
\begin{align*}
U_{1,h}^\varepsilon(b_\varepsilon)
&
=\Big|\int_0^t\(P^{Y^\varepsilon}_s(b_\varepsilon\cdot\nabla u)
(s,\cdot,\mu^\varepsilon_s)(x)-P^{Y^\varepsilon}_{\pi_h(s)}
(b_\varepsilon\cdot\nabla u)(s,\cdot,\mu^\varepsilon_s)(x)\){\mathord{{\rm d}}} s\Big|\\
&
\lesssim \int_0^t \[(h^{\alpha}(\pi_h(s))^{-\alpha})\wedge1\]
(\pi_h(s))^{-\frac{d}{2p_0}}|\mspace{-3mu}|\mspace{-3mu}| (b_\varepsilon\cdot\nabla u)(s,\cdot,\mu^\varepsilon_s)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}{\mathord{{\rm d}}} s\\
&
\lesssim h^\alpha\int_0^t (\pi_h(s))^{-\alpha}
(\pi_h(s))^{-\frac{d}{2p_0}}(t-s)^{-\frac12}{\mathord{{\rm d}}} s,
\end{align*}
where $\alpha=1/2-d/(2p_0)$.
Noting that $\pi_h(s)=s$ for $s\leqslant h$, $\pi_h(s)\leqslant s$ for all $s\in[0,T]$ and $1-\alpha-d/(2p_0)-1/2=0$, one sees that
\begin{align*}
U_{1,h}^\varepsilon(b_\varepsilon)\lesssim h^{\frac12-\frac{d}{2p_0}}.
\end{align*}
For $U_{2,h}^\varepsilon(b_\varepsilon)$, by \eqref{BB06}, \eqref{NN04} and {\bf (H$^1_b$)}, we have
\begin{align*}
U_{2,h}^\varepsilon(b_\varepsilon)\lesssim \int_0^t (\pi_h(s))^{-\frac{d}{2p_0}}\|\mu^\varepsilon_s-\mu^\varepsilon_{\pi_h(s)}\|_{var}(t-s)^{-\frac12}{\mathord{{\rm d}}} s.
\end{align*}
It follows from \eqref{LL02} with $\delta=\frac{d}{2p_0}$ that
\begin{align*}
U_{2,h}^\varepsilon(b_\varepsilon)\lesssim h^{\frac12-\frac{d}{2p_0}}\int_0^t (\pi_h(s))^{-\frac{d}{2p_0}}(t-\pi_h(s))^{-\frac12}{\mathord{{\rm d}}} s\lesssim h^{\frac12-\frac{d}{2p_0}},
\end{align*}
provided by $p_0>d$.
Finally, in view of \eqref{BB06} and \eqref{FStar00}, by $p_0<\infty$, we have
\begin{align*}
U_{3,h}^\varepsilon(b_\varepsilon)\lesssim h^{\frac12-\frac{d}{2p_0}}\int_0^t (\pi_h(s))^{-\frac{d}{2p_0}}(t-\pi_h(s))^{-1+\frac{d}{2p_0}}{\mathord{{\rm d}}} s
\lesssim h^{\frac12-\frac{d}{2p_0}}
\end{align*}
and obtain that
\begin{align}\label{NN06}
\sum_{i=1}^3\(U_{i,h}^\varepsilon(b_\varepsilon)+U_{i,h}^\varepsilon(\bar{b})\)\lesssim h^{\frac12-\frac{d}{2p_0}}.
\end{align}

{\bf(Step 2)} First of all, let $M:=[t/h]$. Without loss of generality we assume that $M=t/h\in{\mathbb N}$ and notice that
\begin{align*}
{\mathcal E}^\varepsilon_h\leqslant&\Big|{\mathbb E}\int_{0}^{h}\Big(b_\varepsilon(s,Y^\varepsilon_{s}(x),\mu^\varepsilon_{s})-b_0(Y^\varepsilon_{s}(x),\mu^\varepsilon_{s})\Big)\cdot\nabla u(s,Y^\varepsilon_{s}(x)){\mathord{{\rm d}}} s\Big|\\
&
+ \Big|\sum_{k=1}^{M-1}{\mathbb E}\int_{kh}^{(k+1)h}
\Big(b_\varepsilon(s,Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})
-b_0(Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})\Big)\cdot\nabla u(kh,Y^\varepsilon_{kh}(x)){\mathord{{\rm d}}} s\Big|\\
:=&{\mathcal E}_1+{\mathcal E}_2.
\end{align*}
From \eqref{BB06} and \eqref{NN04},
\begin{align*}
{\mathcal E}_1\lesssim \int_{0}^{h} s^{-\frac{d}{2p_0}}(t-s)^{-\frac12}{\mathord{{\rm d}}} s\lesssim h^{\frac12-\frac{d}{2p_0}}.
\end{align*}
By \eqref{NN04} and a change of variable, one sees that
\begin{align*}
{\mathcal E}_2&\lesssim \sum_{k=1}^{M-1}(t-kh)^{-\frac12}\Big|{\mathbb E}\int_{kh}^{(k+1)h}\Big(b_\varepsilon(s,Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})-b_0(Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})\Big){\mathord{{\rm d}}} s\Big|\\
&\lesssim \sum_{k=1}^{M-1}(t-kh)^{-\frac12}\Big|\varepsilon{\mathbb E}\int_{kh/\varepsilon}^{(k+1)h/\varepsilon}\Big(b(s,Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})-\bar{b}(Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})\Big){\mathord{{\rm d}}} s\Big|.
\end{align*}
Based on the assumption \eqref{in:00} and \eqref{BB06}, we have
\begin{align*}
{\mathcal E}_2
&
\lesssim h\sum_{k=1}^{M-1}(t-kh)^{-\frac12}\omega\left(\frac{h}{\varepsilon}\right){\mathbb E} H(Y^\varepsilon_{kh}(x),\mu^\varepsilon_{kh})\\
&
\lesssim h\sum_{k=1}^{M-1}(t-kh)^{-\frac12}\omega\left(\frac{h}{\varepsilon}\right) (kh)^{-\frac{d}{2p_0}} \sup_\mu|\mspace{-3mu}|\mspace{-3mu}| H(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\\
&
\lesssim \omega\left(\frac{h}{\varepsilon}\right) \int_h^t (t-\pi_h(s))^{-\frac12}(\pi_h(s))^{-\frac{d}{2p_0}}{\mathord{{\rm d}}} s
\lesssim \omega\left(\frac{h}{\varepsilon}\right)
\end{align*}
and obtain that
\begin{align*}
|{\mathbb E} \varphi(X^\varepsilon_t)-{\mathbb E} \varphi(X_t)|\lesssim\(h^{\frac12-\frac{d}{2p_0}}+\omega\left(\frac{h}{\varepsilon}\right)\)+\int_0^ts^{-\frac{d}{2p_0}}(t-s)^{-\frac12}\|\mu^\varepsilon_s-\mu_s\|_{var}{\mathord{{\rm d}}} s
\end{align*}
provided by \eqref{NN03}, \eqref{NN05} and \eqref{NN06}.
Finally, taking the supremum of al $\varphi\in C^\infty_b$ with $\|\varphi\|_\infty=1$ and
by Gronwall's inequality of Volterra type (see \cite[Example 2.4]{Zhang2010}), we complete the proof.

\end{proof}

\section{Proof of Theoerm \ref{in:Main} and \ref{in:Main2}}\label{Section5}
In this section, we consider the process $(X^\varepsilon,W,X)$ on probability space $(\Omega,{\mathscr F},({\mathscr F}_t)_{t\ge0}, {\mathbb P})$ which satisfies the following systems in ${\mathbb R}^d$:
\begin{align*}
X^\varepsilon_t=\xi+\int_0^t b(\frac{s}{\varepsilon},X^\varepsilon_s,\mu^\varepsilon_s){\mathord{{\rm d}}} s+\int_0^t\sigma(X^\varepsilon_s){\mathord{{\rm d}}} W_s
\end{align*}
and
\begin{align*}
X_t=\xi+\int_0^t \bar{b}(X_s,\mu_s){\mathord{{\rm d}}} s+\int_0^t\sigma(X_s){\mathord{{\rm d}}} W_s,
\end{align*}
where $W$ is a standard $d$-dimensitional Brownian Motion, $\mu^\varepsilon_t$ and $\mu_t$ are the distributions of $X^\varepsilon_t$ and $X_t$ respectively and $(b,\bar{b},\sigma)$ satisfies the conditions {\bf (H$_b^1$)}-{\bf (H$_b^2$)} and {\bf (H$_\sigma$) }. For the unity of the notations, we set
\begin{align*}
X^0:=X,\ \ b_\varepsilon(t):=b(t/\varepsilon),\ \ \text{and} \ \ b_0:=\bar{b}.
\end{align*}

\begin{proof}[Proof of Theoerm \ref{in:Main}]
Set
\begin{align*}
B(t,x):=b_0(x,\mu_t)
\end{align*}
and consider the following backward parabolic PDE
$$
\partial_tu+a_{ij}\partial_i\partial_j u-\lambda u+B\cdot\nabla u+B=0, ~~t\in[0,T],\ \ u(T)=0.
$$
Noting that $|\mspace{-3mu}|\mspace{-3mu}| B|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}
\leqslant\sup\limits_\mu|\mspace{-3mu}|\mspace{-3mu}| \bar{b}(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}<\infty $, by Lemma \ref{lem210}, taking $\lambda$ large enough, there is a unique solution $u$ in the sense of Definition \ref{Pre4:Def} satisfying
\begin{align*
\|\nabla u\|_{{\mathbb L}^\infty_T}\leqslant \frac{1}{2}
\end{align*}
and for any $2/q+d/p_0<1$,
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}|\nabla^2 u|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}_q(T)}\leqslant C,
\end{align*}
  which implies that for any $\lambda>0$
  \begin{align}\label{Khau1}
  \sup_{\varepsilon\ge0}{\mathbb E} \exp\left(\lambda \int_0^T|\nabla^2 u(t,X^{\varepsilon}_t)|^2{\mathord{{\rm d}}} t\right)<\infty,
\end{align}
where $X^0:=X$,
  provided by \eqref{Kh02}. Moreover, by \eqref{HH02}, for all $s,t\in[0,T]$,
  \begin{align}\label{GG051}
  \|\nabla u(t)-\nabla u(s)\|_\infty \lesssim |t-s|^{1/2-d/(2p_0)}.
\end{align}
Define
$$
\Phi_t(x):=x+u(t,x)
$$
and
$$
Y^{\varepsilon}_t:=\Phi_t({X}^{\varepsilon}_t),\ \ Y_t:=\Phi_t(X_t).
$$
Then $\Phi_t$ is a $C^1$-diffeomorphism (see Remark \ref{Re002}) for any $t\in[0,T]$ with
\begin{align}\label{ZZ031}
\|\nabla\Phi\|_{{\mathbb L}^\infty_T}+\|\nabla \Phi^{-1}\|_{{\mathbb L}^\infty_T}\leqslant 4.
\end{align}
By generalized It\^o's formula \eqref{GIF001}, we have
\begin{align*}
{\mathord{{\rm d}}} Y_t=\lambda u(t,{X}_t){\mathord{{\rm d}}} t+(\sigma^*\nabla\Phi_t)({X}_t){\mathord{{\rm d}}} W_t
\end{align*}
and
\begin{align*}
{\mathord{{\rm d}}} Y^{\varepsilon}_t=\lambda u(t,X^\varepsilon_t){\mathord{{\rm d}}} t+\(b_\varepsilon(t,X^\varepsilon_t,\mu^\varepsilon_t)-b_0(X^\varepsilon_t,\mu_t)\)\cdot\nabla \Phi_t\(X^{\varepsilon}_t\){\mathord{{\rm d}}} t+(\sigma^*\nabla\Phi_t(X^\varepsilon_t)){\mathord{{\rm d}}} W_t,
\end{align*}
where $\sigma^*$ is the transpose of $\sigma$.
It follows from \eqref{ZZ031} that for any $t\in[0,T]$,
\begin{align*}
|X^\varepsilon_t-X_t|^2\lesssim |Y^\varepsilon_t-Y_t|^2
&
\lesssim \lambda\|\nabla u\|_{{\mathbb L}^\infty_T}^2\int_0^t |X^\varepsilon_s-X_s|^2{\mathord{{\rm d}}} s\\
&
+\left[\int_0^t\((\sigma^*\nabla\Phi_s)(X^\varepsilon_s)-(\sigma^*\nabla\Phi_s)(X_s)\){\mathord{{\rm d}}} W_s\right]^2\\
&
+\Big|\int_0^t\(b_\varepsilon(s,X^\varepsilon_s,\mu^\varepsilon_s)-b_0(X^\varepsilon_s,\mu_s)\)\cdot\nabla \Phi_s\(X^{\varepsilon}_s\){\mathord{{\rm d}}} s\Big|^2.
\end{align*}
Set
\begin{align*}
A^{\varepsilon}_t:=&\int_0^t\left({\mathcal M}(\nabla^2 u)(s,X_s)+{\mathcal M}(\nabla^2 u)(s,X^\varepsilon_s)+\|\nabla u\|_{{\mathbb L}^\infty_T}\right)^2{\mathord{{\rm d}}} s\\
&+\int_0^t\left({\mathcal M}(\nabla \sigma)(X_s)+{\mathcal M}(\nabla \sigma)(X^\varepsilon_s)+\|\sigma\|_\infty\right)^2{\mathord{{\rm d}}} s
\end{align*}
and
\begin{align*}
\eta^\varepsilon_t:=\Big|\int_0^t\(b_\varepsilon(s,X^\varepsilon_s,\mu^\varepsilon_s)-b_0(X^\varepsilon_s,\mu_s)\)\cdot\nabla \Phi_s\(X^{\varepsilon}_s\){\mathord{{\rm d}}} s\Big|^2.
\end{align*}
Then, by \eqref{Pre1:Mix}, \eqref{Khau1} and {\bf (H$_\sigma$)}, we have
\begin{align}\label{Coffee}
\sup_\varepsilon{\mathbb E} \exp(A^{\varepsilon}_T)<\infty.
\end{align}
Notice that by \eqref{Pre1:Mnf}
\begin{align*}
&\quad
\left[\int_0^t\((\sigma^*\nabla\Phi)(X^\varepsilon_s)-(\sigma^*\nabla\Phi)(X_s)\){\mathord{{\rm d}}} W_s\right]^2\\
&
\leqslant \int_0^t|X^\varepsilon_s-X_s|^2{\mathord{{\rm d}}} A^{\varepsilon}_s+M^\varepsilon_t,
\end{align*}
where $M^\varepsilon$ is a martingale.
Summary, we have
\begin{align*}
|X^\varepsilon_t-X_t|^2\lesssim\int_0^t|X^\varepsilon_s-X_s|^2{\mathord{{\rm d}}} s+\int_0^t|X^\varepsilon_s-X_s|^2{\mathord{{\rm d}}} A^{\varepsilon}_s+M^\varepsilon_t+\eta_t^\varepsilon.
\end{align*}
Hence,
 by \eqref{Coffee} and Stochastic Gronwall's inequality (see Lemma 2.8 in \cite{Zh-Zh1} or \cite{Sc13}), one sees that for any $\ell\in(0,1)$,
\begin{align*}
{\mathbb E}\(\sup_{t\in[0,T]}|X^\varepsilon_t-X_t|^{2\ell}\)\lesssim \({\mathbb E}[\sup_{t\in[0,T]} \eta_t^\varepsilon]\)^\ell.
\end{align*}

Combining \eqref{GG051}, \eqref{ZZ031} and \eqref{EE01}, we have
\begin{align}\label{GG07}
{\mathbb E}[\sup_{t\in[0,T]}\eta_t^\varepsilon]\lesssim
&
{\mathbb E}\int_0^T\Big|b_0(X^\varepsilon_s,\mu^\varepsilon_s)-b_0(X^\varepsilon_s,\mu_s)\Big|^2{\mathord{{\rm d}}} s+ \inf_{h>0}\Big(h^{1-\delta}+h^{1-\frac{d}{p_0}}+\omega\left(\frac{h}{\varepsilon}\right)\Big).
\end{align}
Taking $\delta<d/p_0$ in \eqref{GG07}, from \eqref{Kry02} and \eqref{LL00},
for any $2/q+d/p_0<1$, one sees that
\begin{align*}
{\mathbb E}[\sup_{t\in[0,T]}\eta_t^\varepsilon]
&
\lesssim\Big(\int_0^T\|\mu^\varepsilon_s-\mu_s\|_{var}^q{\mathord{{\rm d}}} s\Big)^{2/q} +\inf_{h>0}\Big(h^{1-\frac{d}{p_0}}
+\left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2\Big)\\
&
\lesssim\inf_{h>0}\Big(h^{1-\frac{d}{p_0}}
+\left(\omega\left(\frac{h}{\varepsilon}\right)\right)^2\Big)
\end{align*}
and this completes the proof.

\end{proof}

In the rest of this section, we assume that
\begin{align*}
b_\varepsilon(t,x,\mu)=b_\varepsilon(t,x),\quad b_0(x,\mu)=b_0(x)
\end{align*}
and prove Theorem \ref{in:Main2}. The method is totally the same as the one of Theorem \ref{in:Main}, except for using the elliptic equation to construct the Zvonkin's transformation.

\begin{proof}[Proof of Theorem \ref{in:Main2}]
Consider the following elliptic PDE
\begin{align}\label{Ell00}
a_{ij}\partial_i\partial_j u-\lambda u+b_0\cdot\nabla u+b_0=0.
\end{align}
Noting that $|\mspace{-3mu}|\mspace{-3mu}| b_0|\mspace{-3mu}|\mspace{-3mu}|_{\widetilde{{\mathbb L}}^{p_0}(T)}\leqslant\|b\|_{L^\infty({\mathbb R}_+;\widetilde{L}^{p_0})} $, by \eqref{CC02}, taking $\lambda$ large enough, we have
\begin{align*
\|\nabla u\|_{\infty}\leqslant \frac{1}{2}
\end{align*}
and
\begin{align*}
|\mspace{-3mu}|\mspace{-3mu}|\nabla^2 u|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\leqslant C \|b_0\|_{\widetilde{L}^{p_0}},\quad \forall p_0>d
\end{align*}
 It follows that for any $\lambda>0$
  \begin{align}\label{Khau}
  \sup_{\varepsilon\ge0}{\mathbb E} \exp\left(\lambda \int_0^T|\nabla^2 u(X^{\varepsilon}_t)|^2{\mathord{{\rm d}}} t\right)<\infty
\end{align}
by \eqref{Kh02}.
Define
$$
\Phi(x):=x+u(x)
$$
and
$$
Y^{\varepsilon}_t:=\Phi({X}^{\varepsilon}_t),\ \ Y_t:=\Phi(X_t).
$$
Then $\Phi$ is a $C^1$-diffeomorphism
By generalized It\^o's formula \eqref{GIF001} again, we have
\begin{align*}
{\mathord{{\rm d}}} Y_t=\lambda u({X}_t){\mathord{{\rm d}}} t+(\sigma^*\nabla\Phi)({X}_t){\mathord{{\rm d}}} W_t
\end{align*}
and
\begin{align*}
{\mathord{{\rm d}}} Y^{\varepsilon}_t=\lambda u(X^\varepsilon_t){\mathord{{\rm d}}} t+\(b_\varepsilon(t,X^\varepsilon_t)-b_0(X^\varepsilon_t)\)\cdot\nabla \Phi\(X^{\varepsilon}_t\){\mathord{{\rm d}}} t+(\sigma^*\nabla \Phi)(X^\varepsilon_t){\mathord{{\rm d}}} W_t.
\end{align*}
Then, we have
\begin{align*}
|X^\varepsilon_t-X_t|^2\lesssim \int_0^t |X^\varepsilon_s-X_s|^2{\mathord{{\rm d}}} A^\varepsilon_s+M^\varepsilon_t+\eta_t^\varepsilon,
\end{align*}
where $(M^\varepsilon_t)_{t\ge0}$ is a martingale,
\begin{align*}
A^{\varepsilon}_t=&t+\int_0^t\left({\mathcal M}(\nabla^2 u)(X_s)+{\mathcal M}(\nabla^2 u)(X^\varepsilon_s)+\|\sigma\|_\infty\right)^2{\mathord{{\rm d}}} s\\
&+\int_0^t\left({\mathcal M}(\nabla \sigma)(X_s)+{\mathcal M}(\nabla \sigma)(X^\varepsilon_s)+\|\nabla u\|_{\infty}\right)^2{\mathord{{\rm d}}} s
\end{align*}
and
\begin{align*}
\eta^\varepsilon_t=\Big|\int_0^t\(b(s/\varepsilon,X^\varepsilon_s)-\bar {b}(X^\varepsilon_s)\)\cdot\nabla \Phi\(X^{\varepsilon}_s\){\mathord{{\rm d}}} s\Big|^2.
\end{align*}
Then, in view of \eqref{Pre1:Mix} and \eqref{Khau}, we have
\begin{align*}
\sup_\varepsilon{\mathbb E} \exp(A^{\varepsilon}_T)<\infty,
\end{align*}
which implies that for any $\ell\in(0,1)$,
\begin{align*}
{\mathbb E}\(\sup_{t\in[0,T]}|X^\varepsilon_t-X_t|^{2\ell}\)\lesssim \({\mathbb E}\sup_{t\in[0,T]} \eta_t^\varepsilon\)^\ell
\end{align*}
provided by Stochastic Gronwall's inequality and we complete the proof by \eqref{EE01} with $\alpha=1$.
\end{proof}

\section*{}
\subsection*{Acknowledgments}
The first author would like to acknowledge the
warm hospitality of Bielefeld University. And
we would also like to thank Dr. Chengcheng Ling for many useful discussions.

\begin{appendix}
\renewcommand{\thetable}{A\arabic{table}}
\numberwithin{equation}{section}

\section{}
In appendix, we check the coefficients in Example \ref{Ex2} satisfy our assumptions {\bf (H$_b^1$)} and {\bf (H$_b^2$)}. Since {\bf (H$_b^1$)} is obviously,
we only check the assumption {\bf (H$_b^2$)} by showing the following result.
\begin{lemma}\label{lem61}
Let
$$
b(t):=\int_{\R} B(\sin(\xi t))\nu({\mathord{{\rm d}}} \xi).
$$
and
$$
\bar{b}:=\Big(\frac{1}{2\pi}\int_0^{2\pi}
B(\sin(\xi \tau)){\mathord{{\rm d}}}\tau\Big)\nu(\R\setminus\{0\})
+B(0)\nu(\{0\}),
$$
where $B: [-1,1]$ is measurable and $\nu$ is a finite measure on ${\mathbb R}$.
Assume that
\begin{align*}
|B(u)|\leqslant C_B,\quad \forall u\in[-1,1].
\end{align*}
Then, for any $t,T\in{\mathbb R}_+$,
\begin{align*}
\Big|\frac{1}{T}\int_{t}^{t+T}\(b(s)-\bar{b}\){\mathord{{\rm d}}} s\Big|\leqslant \frac{4\pi C_F}{T}\int_{{\mathbb R}\setminus\{0\}}\frac{\nu({\mathord{{\rm d}}} \xi)}{|\xi|}.
\end{align*}

\end{lemma}
\begin{remark}\rm
Notice that in Example \ref{Ex2}
\begin{align*}
|F(u,x)|\leqslant|F(t,0)|+|F(t,x)-F(t,x)|\leqslant L_F+L_F|x|
\end{align*}
provided by \eqref{in:EX2}, which implies that
\begin{align*}
|F\(u,\int_{{\mathbb R}^d}\phi(x,y)\mu({\mathord{{\rm d}}} y)\)|\leqslant L_F(1+\int_{{\mathbb R}^d}|\phi(x,y)|\mu({\mathord{{\rm d}}} y)).
\end{align*}
Hence, by Lemma \ref{lem61}, we see that coefficients in Example \ref{Ex2} satisfy {\bf (H$_b^2$)} with
$$
\omega(T)=\frac{4\pi L_F}{T}\int_{{\mathbb R}\setminus\{0\}}\frac{\mu({\mathord{{\rm d}}} \xi)}{|\xi|}
$$
and $H(x,\mu)=1+\int_{{\mathbb R}^d}|\phi(x,y)|\mu({\mathord{{\rm d}}} y)$, where
\begin{align*}
\sup_\mu|\mspace{-3mu}|\mspace{-3mu}| H(\cdot,\mu)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\leqslant \|1\|_\infty+\int_{{\mathbb R}^d}|\mspace{-3mu}|\mspace{-3mu}|\phi(\cdot,y)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}\mu({\mathord{{\rm d}}} y)\le1+\sup_{y}|\mspace{-3mu}|\mspace{-3mu}|\phi(\cdot,y)|\mspace{-3mu}|\mspace{-3mu}|_{p_0}.
\end{align*}

\end{remark}
\begin{proof}[Proof of Lemma \ref{lem61}]
When $\int_{{\mathbb R}\setminus\{0\}}\frac{\nu({\mathord{{\rm d}}} \xi)}{|\xi|}=\infty$, it is trivial. So we assume that $\int_{{\mathbb R}\setminus\{0\}}\frac{\nu({\mathord{{\rm d}}} \xi)}{|\xi|}<\infty$. Firstly, one sees that
\begin{align*}
{\mathscr I}
&
:=\left|\frac{1}{T}\int_{t}^{t+T}\(b(s)
-\bar{b}\){\mathord{{\rm d}}} s\right|\\
&
=\left|\frac{1}{T}\int_{t}^{t+T}\int_{{\mathbb R}\setminus\{0\}} B(\sin(\xi s))\nu({\mathord{{\rm d}}}\xi){\mathord{{\rm d}}} s-\frac{1}{2\pi}\int_t^{t+2\pi} B(\sin(\tau))\nu({\mathbb R}\setminus\{0\}){\mathord{{\rm d}}} \tau\right|\\
&=\left|\int_{{\mathbb R}\setminus\{0\}} \Big[\frac{1}{T}\int_{t}^{t+T}B(\sin(\xi s)){\mathord{{\rm d}}} s-\frac{1}{2\pi}\int_t^{t+2\pi} B(\sin(\tau)){\mathord{{\rm d}}} \tau\Big]\nu({\mathord{{\rm d}}}\xi)\right|,
\end{align*}
provided by Fubini's theorem. From a change of variable, we have
\begin{align*}
{\mathscr I}
=&\left|\int_{{\mathbb R}\setminus\{0\}} \Big[\frac{1}{T|\xi|}\int_{t\xi}^{(t+T)\xi}B(\sin s){\mathord{{\rm d}}} s-\frac{1}{2\pi}\int_{t}^{t+2\pi} B(\sin \tau){\mathord{{\rm d}}} \tau\Big]\nu({\mathord{{\rm d}}}\xi)\right|,
\end{align*}
where $\int_a^b:=-\int_b^a$ if $a>b$.
Denote by
\begin{align*}
G:=\int_{0}^{2\pi} B(\sin s){\mathord{{\rm d}}} s=\int_{t}^{t+2\pi} B(\sin s){\mathord{{\rm d}}} s,\quad\forall t\in{\mathbb R}.
\end{align*}
Then,  noting that $s\to\sin s$ has a period $2\pi$, we have
\begin{align*}
\int_{t\xi}^{(t+T)\xi}B(\sin s){\mathord{{\rm d}}} s&=\Big[\frac{T|\xi|}{2\pi}\Big]G+\int_{\mbox{\rm sgn}(\xi)\left[\frac{T|\xi|}{2\pi}\right]2\pi+t\xi}^{T\xi+t\xi} B(\sin s){\mathord{{\rm d}}} s\\
&:=\Big[\frac{T|\xi|}{2\pi}\Big]G+H_t(\xi)
\end{align*}
where $\mbox{\rm sgn}(\xi):=\xi/|\xi|$, which implies that
\begin{align*}
{\mathscr I}
=&\left|\int_{{\mathbb R}\setminus\{0\}} \Big[\frac{1}{T|\xi|}\Big(\Big[\frac{T|\xi|}{2\pi}\Big]G+H_t(\xi)\Big)-\frac{1}{2\pi}G\Big]\nu({\mathord{{\rm d}}}\xi)\right|\\
\leqslant&\int_{{\mathbb R}\setminus\{0\}} \Big|\frac{1}{T|\xi|}\Big[\frac{T|\xi|}{2\pi}\Big]-\frac{1}{2\pi}\Big|\nu({\mathord{{\rm d}}}\xi)G+\int_{{\mathbb R}\setminus\{0\}} \frac{1}{T|\xi|}H_t(\xi)\nu({\mathord{{\rm d}}}\xi).
\end{align*}
Notice that
\begin{align*}
\Big|\frac{1}{T|\xi|}\Big[\frac{T|\xi|}{2\pi}\Big]-\frac{1}{2\pi}\Big|=\frac{1}{T|\xi|}\Big|\Big[\frac{T|\xi|}{2\pi}\Big]-\frac{T|\xi|}{2\pi}\Big|\leqslant \frac{1}{T|\xi|}
\end{align*}
and
\begin{align*}
G\vee H_t(\xi)\leqslant \int_0^{2\pi}|B(\sin s)|{\mathord{{\rm d}}} s\leqslant 2\pi C_B.
\end{align*}
Therefore, we have
\begin{align*}
{\mathscr I}\leqslant \frac{4\pi C_B}{T}\int_{{\mathbb R}\setminus\{0\}}\frac{\nu({\mathord{{\rm d}}} \xi)}{|\xi|}
\end{align*}
and complete the proof.

\end{proof}

\end{appendix}

redpajama_set_name: RedPajamaArXiv