\section{Introduction}

Motivated by analyzing modern machine learning algorithms, we provide a strong uniform consistency result with the convergence rate for kernel density estimation (KDE) with general kernels on Riemannian manifolds.
Our main results can be summarized as follows. Let $M$ be a smooth $d$ dimensional compact manifold without boundary isometrically embedded in $\mathbb{R}^p$ through $\iota$. Let $x_1, \cdots, x_n$ be a sequence of points independently and identically (i.i.d) sampled from $M$ following the probability density function $P$ on $M$. 
\begin{enumerate}
\item
Suppose $K(t):\mathbb{R}_{\geq 0} \rightarrow \mathbb{R}$ (not necessary non-negative) is a bounded {\em Riemann integrable} kernel function with a proper decay rate.  Suppose $\int_{ \mathbb{R}^d} K(\|\textbf{v}\|_{\mathbb{R}^d})d\textbf{v}=1$. Consider the kernel density estimator at $x\in M$ as
\begin{align}\label{Definition Kn isotropic}
K_n(x)=\frac{1}{n \epsilon_n^d}\sum_{i=1}^n K\left(\frac{\|\iota(x_i)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon_n}\right)\,.
\end{align}
Under a suitable relationship between $\epsilon_n$ and $n$, we show that $\sup_{x \in M}|K_n(x)-\mathbb{E} K_n(x)| \rightarrow 0$ a.s. and we will provide the convergence rate. If we further assume that the probability density function $P$ is H\"older continuous on $M$, then we have $\sup_{x \in M}|\mathbb{E} K_n(x)-P(x)| \rightarrow 0$ with a convergence rate.  {See Chapter \ref{Section:KDE2}. This setup is a direct generalization of most traditional setups by considering the manifold and possibly irregular function as an isotropic kernel.}

\item{
Suppose $K(\textbf{v}):\mathbb{R}^d \rightarrow \mathbb{R}$ (not necessary non-negative) is a bounded {\em Riemann integrable} kernel function with compact support.  Suppose $\int_{ \mathbb{R}^d} K(\textbf{v})d\textbf{v}=1$. For any $x$, we have a {diffeomorphism (or chart)} $\Phi_x$ of $\iota(M)$ such that $\Phi_x(\textbf{v}): B^{\mathbb{R}^d}_r(0) \rightarrow \iota(M)$ with $\Phi_x(0)=\iota(x)$, where $B^{\mathbb{R}^d}_r(0)$ is the open Euclidean ball of radius $r$ centered at $0$. Let $U_x(0)$ be the volume form of $\iota(M)$ in the chart $\Phi_x$ at $\textbf{v}=0$. For $x \in M$, we define the kernel density estimator at $x$ to be
\begin{align}
K_{n}(x)=\frac{1}{n \epsilon_n^d U_x(0)}\sum_{i=1}^n K_{\epsilon_n}(x,x_i)\,,
\end{align}
where  $K_{\epsilon_n}(x,x_i)=K(\frac{\Phi^{-1}_x(\iota(x_i))}{\epsilon_n})$ if $\iota(x_i) \in \Phi_x(B^{\mathbb{R}^d}_r(0))$ and $K_{\epsilon_n}(x,x_i)=0$ otherwise. We assume that for all $x$, we have a uniform control of the bi-Lipschitz constant of $\Phi_x$ and the $C^1$ norm of the volume form of $\iota(M)$ in the chart $\Phi_x$ at $\textbf{v}=0$. We further assume that the probability density function $P$ is H\"older continuous on $M$. Under a suitable relationship between $\epsilon_n$ and $n$, we show that $\sup_{x \in M}|\mathbb{E} K_n(x)-P(x)| \rightarrow 0$ with a convergence rate.  
}{See Chapter \ref{Section:KDE3}. This setup is less considered in the traditional setup, but an emerging one inspired by analyzing modern manifold learning algorithms. Here, the kernel $K_{\epsilon}(x,y)$ is in general not translational invariant and irregular.}

\item Suppose $K_{\epsilon}: \iota(M) \times \iota(M) \rightarrow \mathbb{R}_{\geq 0}$ is a sequence of bounded {\em Lebesgue integrable} kernel with compact support indexed by $\epsilon>0$.  Moreover, $\int_{M} K_{\epsilon}(\iota(x),\iota(y))dV(y)=1$, where $dV$ is the volume form on $M$. Let
\begin{align}
K_{n,\epsilon}(x)=\frac{1}{n}\sum_{i=1}^n K_{\epsilon}(\iota(x),\iota(x_i))\,.
\end{align}
Under a suitable relationship between $\epsilon$ and $n$, we show that $K_{n,\epsilon}(x) \rightarrow P(x)$ a.s. in the $L^1$ sense under different conditions on $P$. {See Chapter \ref{Section:KDE4}. This is the most general setup, where the kernel is in general not translational invariant and only Lebesgue integrable, and only $L^1$ convergence is shown.}
\end{enumerate}

KDE on manifolds shows up{, sometimes obviously and sometimes occultly, in manifold learning algorithms. Take the} nonlinear dimension reduction (NDR) {as an example}.
The purpose of an NDR algorithm is to reconstruct the dataset on a manifold of a high dimensional space in some low dimensional space, while {preserving} the geometric or topological structure of the dataset.  However, the {non-uniform} distribution of the dataset on the manifold is one of the major obstacles to the exploration of the geometric or topological structure of the dataset. Hence, {it is intuitive that} KDE plays an essential role to eliminate the impact of the {non-uniform} distribution of the dataset when one develops an NDR algorithm. 

{There are two main ways KDE might appear in a manifold learning algorithm.
The first one comes from a direct design aiming to eliminating the impact of non-uniform distribution.}
For instance, the $\alpha=1$ renormalization in the diffusion map \cite{coifman2006diffusion} is the step to estimate the density of the dataset based on the kernel chosen by the user. {For KDE of this kind of appearance, since the manifold is usually unknown, the kernel is usually designed in the ambient Euclidean space, and we need to handle the intricate interaction between the intrinsic and extrinsic geometries when we analyze its behavior. Usually, researchers consider isotopic (that is, the kernel is defined as $K(x,y):=f(\|x-y\|)$ for a ``proper'' function $f$ defined on $\mathbb{R}_{\geq 0}$ for all $x,y\in \iota(M)$), or sometimes anisotropic (that is, the kernel is defined as $K(x,y):=k(x-y)$ for a ``proper'' function $k$ defined on $\mathbb{R}^p$ for all $x,y\in \iota(M)$), kernels that are symmetric with some nice properties, but rarely consider consider non-translational invariant (that is, the kernel $K$ is defined on $\mathbb{R}^p\times \mathbb{R}^p$ and cannot be defined on $x-y$ for all $x,y\in \mathbb{R}^p$) kernels with undesired properties like irregularity. We will discuss the isotropic case in Section \ref{Section:KDE2}.}


{The second one is not obvious from directly reading the algorithms, and the involvement of KDE was probably not the focus when designing the algorithm. However, it appears naturally when we analyze these algorithms under the manifold setup. The associated kernel of this kind is not designed by the user but appears intrinsically when the algorithm is applied, is usually defined on the manifold, and is probably not translational invariant, dependent on $n$ and other quantities, and irregular}. The widely applied locally linear embedding (LLE) algorithm \cite{Roweis_Saul:2003} is {a good} example. It is shown in \cite[Corollary 3.1]{wu2018think} that under the manifold setup, the LLE kernel, {and hence KDE,} is adaptive to the geometry of the dataset and depends on the regularization chosen in the barycentric coordinate evaluation. Therefore, the kernel is {non-translational invariant} and irregular. {Moreover,} the kernel's behavior near the boundary of a manifold is dramatically different from that away from the boundary \cite{wu2018locally}. As shown in \cite{wu2018think} and \cite{wu2018locally}, KDE is a crucial step towards the asymptotic analysis of LLE. {Getting a deeper theoretical justification of LLE relies on studying the spectral behavior of the integral operator corresponding to the LLE kernel, where a $L^\infty$ sense convergence of KDE with a general kernel is involved. {We will discuss this case in Sections \ref{Section:KDE3} and \ref{Section:KDE4}.}


Therefore, although the kernels considered in this work might seem not practical from the traditional perspective, it emerges naturally when we analyze the modern machine learning algorithms, which evolves fast in the past decades. While nowadays LLE might be considered a ``traditional'' NDR algorithm, its non-obvious involvement of KDE and intricate kernel structure suggests the possibility that we might encounter more challenging kernels when we analyze ``modern'' machine learning algorithms.}


Our results could be viewed as a generalization of various existing results.
The traditional KDE on the Euclidean space was first considered by M. Rosenblatt \cite{rosenblatt1956remarks} and E. Parzen \cite{parzen1962estimation}.  Let $P$ be a probability density function on $\mathbb{R}^d$. Let $x_1, \cdots, x_n$ be a sequence of i.i.d samples from $\mathbb{R}^d$ based on $P$. Let $\epsilon_n$ be a sequence of numbers such that $\epsilon:=\epsilon_n \rightarrow 0$ and $n\epsilon^d_n \rightarrow \infty$ as $n \rightarrow \infty$. 
For $x \in \mathbb{R}^d$, the kernel density estimator at $x$ is defined as
\begin{align}
K_{n}(x):=\frac{1}{n \epsilon^d}\sum_{i=1}^n K\left(\frac{x_i-x}{\epsilon}\right)\,, 
\end{align}
where $K(x) \geq 0$ is a bounded real valued function on $\mathbb{R}^d$. After \cite{rosenblatt1956remarks,parzen1962estimation}, a lot of results under various conditions were proposed. We summarize those results that are directly related to our work.

\begin{enumerate}
\item In \cite{devroye1979l1}, the authors consider the KDE on the Euclidean space. It is shown that if $P$ is bounded and $K$ is non-negative and its Lebesgue integration is $1$, then $K_n(x) \rightarrow P(x)$ in $L_1(\mathbb{R}^d)$ a.s. when $n \rightarrow \infty$. It is also shown that if there is no condition on $P$ and $K$ further satisfies certain decay rate, then $K_n(x) \rightarrow P(x)$ in $L_1(\mathbb{R}^d)$ a.s. when $n \rightarrow \infty$.
In \cite{devroye1980strong}, it is shown that if $P$ is uniformly continuous and $K$ is non-negative with certain decay rate and Riemann integration $1$, then $K_n(x) \rightarrow P(x)$ for all $x$, a.s. when $n \rightarrow \infty$. {Density estimation from the Bayesian approach could be found in, for example, \cite{wong2010optional, ma2017adaptive} and many others. }

\item In \cite{gine2002rates} and \cite{einmahl2005uniform}, the authors consider the KDE in the Euclidean space. The authors provide the  {$L^\infty$ convergence rate} from $K_n(x)$ to $\mathbb{E}K_n(x)$ for bounded $P$ and $K$ in the VC class. Compared with the convergence of $K_n(x)$ to $\mathbb{E}K_n(x)$, the regularity of $P$ is necessary for the convergence rate from $\mathbb{E}K_n(x)$ to $P(x)$. The result is applied to the density clustering problem in \cite{rinaldo2010generalized}. The proof of the convergence of $K_n(x)$ to $\mathbb{E}K_n(x)$ relies on the {{\em Rademacher process}}. In this process, it is crucial that the $\epsilon$-covering number of the kernel in the $L^2$ norm is of order $\frac{1}{\epsilon^v}$ for some $v>0$ so that the log of the covering number is an integrable function of $\epsilon$. We mention that the VC class is {\em in some sense} stronger than Riemann integrable in that there exists a bounded  kernel with compact support which is continuous except at one point but not in the VC class. We will discuss this difference more extensively below.


\item There have been several KDE results on closed Riemannian manifold (i.e., compact {\em without} boundary).  In \cite{pelletier2005kernel}, the author works on a closed Riemannian manifold and assumes that $P$ is $C^2$ and $K$ is a non-negative Lebesgue integrable with compact support. Through the Taylor expansion of $P$, the author provides the convergence rate of $\mathbb{E}_P\|K_n-P\|_{L^2(M)}$. In this work, the author uses the intrinsic information, the {\em geodesic distance}, rather than the extrinsic information, the ambient Euclidean distance, to construct the kernel. However, in practice, the geodesic distance information is usually not accessible to researchers. Similarly, the authors in \cite{kim2013geometric} also construct a kernel density estimator by using the intrinsic information of the manifold. Through the exponential map, they pull back the sample points on the manifold to the tangent space and apply KDE for the tangent space (Euclidean space). When the kernel function is regular and the probability density function is $C^4$, they provide an $L^2$ convergence rate result. In contrast, the authors in \cite{ozakin2009submanifold} only use the ambient Euclidean distance to construct a density estimator on the manifold and prove an $L^2$ convergence rate for the {non-negative} $C^1$ kernel functions and $C^2$ probability density function.   In \cite{berenfeld2019density}, the authors study the $L^p$ convergence rate for smooth kernel functions with compact support. The convergence rate  is estimated in terms of the H\"older's exponent of the probability density function and the regularity of the manifold (the H\"older's exponents of the derivatives of the exponential map). In the case when the manifold is compact {\em with} boundary, under a similar assumption as in \cite{pelletier2005kernel}, the authors in \cite{berry2017density} prove an $L^2$ convergence result.   In \cite{chen2019generalized}, the author studies the singular measures on the manifold. The author proves that the KDE is not uniformly consistent but pointwise consistent after rescaling.

\item We mention that the asymptotic convergence of KDE in the {\em pointwise} sense under the manifold setup is widely considered {\em implicitly} in several asymptotic analyses of manifold learning algorithms. In addition to those in the LLE analysis mentioned above \cite{wu2018think,wu2018locally},  see, for example, \cite{belkin2007convergence} and many others. To our knowledge, in most work, since the KDE is not the focus, usually the probability density function and the kernel, if needed, are assumed to behave nicely.  
\end{enumerate}

To our knowledge, our work is the first detailed analysis of the $L^\infty$ convergence rate of KDE with such general conditions on the kernel and the probability density function under the manifold setup. A related but different result under the different kernel assumption is shown in \cite{rinaldo2010generalized}. In that paper, when the dataset is sampled from the manifold, the authors directly generalize [7, Corollary 2.2] to estimate the {\em mollified density function} {that is defined in the ambient space}, where the kernel is in the VC class. The uniform convergence in \cite{rinaldo2010generalized} is in the ambient Euclidean space but not intrinsic to the manifold. As a result, the convergence in \cite{rinaldo2010generalized} might suffer when the ambient space dimension is high. Moreover, it is not clear how to directly ``deconvolve'' the estimated mollified density function to recover the density function on the manifold. Thus, the uniform convergence result under the manifold setup shown in \cite{rinaldo2010generalized} is not comparable with ours. On the other hand, we found that the technical challenge to produce the convergence rate of the KDE on manifolds with kernels classified by their integrability is different from that with kernels in the VC class, and the developed analysis scheme could be useful in other scenarios.

Through the whole paper, we consider the following manifold model. Let $M$ be a smooth $d$ dimensional compact manifold without boundary isometrically embedded in $\mathbb{R}^p$ through $\iota$. Let $P$ be a probability density function on $M$. Let $x_1, \cdots, x_n$ be i.i.d sampled from $M$ based on $P$. We will discuss the KDE with {three} types of kernels: isotropic Riemann integrable kernels, {non-translational invariant Riemann integrable kernels,} and Lebesgue integrable kernels. 
We impose different conditions on the probability density function $P$ and study the convergence behavior of the kernel density estimator.
We use $\textbf{u}=(u_1, \cdots, u_p)$ to denote a vector in $\mathbb{R}^p$ and $\textbf{v}=(v_1, \cdots , v_d)$ to denote a vector in $\mathbb{R}^d$.






\section{Kernel density estimation with isotropic Riemann integrable kernels}\label{Section:KDE2}

We start from stating our assumptions on the kernel.

\begin{Assumption}\label{assumptions on kernel 1}
Suppose $K(t):\mathbb{R}_{\geq 0} \rightarrow \mathbb{R}$ is the kernel function satisfying the following conditions:
\begin{enumerate}[(i)]
\item
$K$ is a bounded function on $\mathbb{R}_{\geq 0}$; that is, $\sup_{t \in \mathbb{R}_{\geq 0}} |K(t)| = K_{\sup}$ for some $K_{\sup}>0$.
\item
$K$ is Riemann integrable on any compact subset $[0,a]$, where $a>0$.
\item
There exists a $\rho>0$, such that if $t \geq \rho$, then $|K(t)| \leq \frac{1}{t^\alpha}$, where $\alpha>d$.
\item
$\int_{ \mathbb{R}^d} K(\|\textbf{v}\|_{\mathbb{R}^d})d\textbf{v}=1$.
\end{enumerate}
\end{Assumption} 

We have a few comments about the above assumptions. First, we do not need the kernel to be non-negative. 
Second, since $K(\|\textbf{u}\|_{\mathbb{R}^p})$ can be regraded as an {\em isotropic} kernel function on the ambient space $\mathbb{R}^p$, we still call the kernel an isotropic one for the density estimation on $M$ despite the fact that it {might not be} isotropic with respect to the intrinsic geometry of the manifold. {The relationship between the intrinsic and extrinsic geometries will be further elaborated in Section \ref{Relationship between intrinsic and extrinsic}.}
Third, note that if $K(t):\mathbb{R}_{\geq 0} \rightarrow \mathbb{R}$ is a bounded function that is Riemann integrable on the compact subsets  $[0,a]$ for any $a$ and $f: A \subset \mathbb{R}^q \rightarrow \mathbb{R}_{\geq 0}$ is a continuous function on a closed rectangle $A$, then it is in general not true that $K \circ f$ is  Riemann integrable on $A$. However, if $f$ is the Euclidean distance function; that is, $f(x)=\|x\|$, we know that $K(\|\textbf{v}\|_{\mathbb{R}^q})$ is Riemann integrable on any closed rectangle $A$ in $\mathbb{R}^q$.
Hence, $\int_{ \mathbb{R}^d} K(\|\textbf{v}\|_{\mathbb{R}^d})d\textbf{v}=1$ in the assumptions can be either understood as an improper Riemann integral or as a Lebesgue integral.

Let $\epsilon_n$ be a sequence of numbers such that $\epsilon_n \rightarrow 0$, as $n \rightarrow \infty$. For $x \in M$, consider the kernel density estimator at $x$ defined in \eqref{Definition Kn isotropic}.
Note that
\begin{align}
\mathbb{E}K_n(x)=\int_{M} \frac{1}{\epsilon_n^d}K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon_n}\right) P(y) dV(y),
\end{align}
where $dV$ is the volume form on $M$. 
We need to emphasize that the integration in the above equation should be understood as {\em Lebesgue integral} on $M$ even though we assume that $K(t)$ is Riemann integrable. In fact, we will see in the later section that under the manifold setup we consider here, for an arbitrary Riemann integrable kernel $K(t)$, $K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon_n}\right)$ {\em may not} be a Riemann integrable function on $M$. 

{In this section,} we expect to show that when $n \rightarrow \infty$, 
\begin{align}
K_n(x) \rightarrow P(x)
\end{align}
a.s. for {\em all} $x$ on $M$ under some convergence rate; that is, $L^\infty$ convergence.
We split this result in two steps, including the variance analysis and the bias analysis. In the bias analysis, under the regularity assumption of $P$, we control $\|\mathbb{E}K_n(x)-P(x)\|_{L^\infty(M)}$ in the deterministic way.
In the variance analysis, we control the variance of $\|\mathbb{E}K_n(x)-K_n(x)\|_{L^\infty(M)}$ in the probabilistic sense, where the regularity assumption of $P$ is not needed. 

\subsection{\textnormal{\textbf{Variance analysis.}}}

We introduce the {\em partition number} for the kernel $K(t)$ satisfying Assumption \ref{assumptions on kernel 1}. 
The partition number plays an essential role in calculating the convergence rate. 

\begin{definition}(Partition number)\label{Definition:Partition number}
Consider the kernel function $K(t)$ satisfying Assumption \ref{assumptions on kernel 1}. Let $D_{lip}$ be a constant that depends only on the $d$-dimensional smooth and compact manifold without boundary embedded in $\mathbb{R}^p$. See \eqref{map J_x} for a precise definition. For any $\gamma>0$, the partition number $N(\gamma)$ is the smallest integer so that the following condition is satisfied. If we partition $[-D_{lip}\gamma^{-\frac{1}{\alpha}}, D_{lip}\gamma^{-\frac{1}{\alpha}}]^d$ uniformly into $N(\gamma)$  cubes $\{\mathcal{Q}_i\}_{i=1}^{N(\gamma)}$, where each $\mathcal{Q}_i$ is of the form $[v_1,v_1+a] \times \cdots \times [v_d,v_d+a]$ for $a>0$, we have 
\begin{align}\label{equation def partition number}
\sum_{i =1}^{N(\gamma)}[M_i (K)-m_i (K)]Vol(\mathcal{Q}_i) <\gamma^2.
\end{align} 
Here, $M_i(K)$ and $m_i(K)$ are  the supremum and infimum of  $K(\|\textbf{v}\|_{\mathbb{R}^d})$ respectively over the cube $\mathcal{Q}_i$. 
\end{definition}

Recall the above-mentioned fact that $K(\|\textbf{v}\|_{\mathbb{R}^d})$ is Riemann integrable over $[-D_{lip}\gamma^{-\frac{1}{\alpha}}, D_{lip}\gamma^{-\frac{1}{\alpha}}]^d$. Hence, the existence of $N(\gamma)$ follows from the Riemann integrability.


\begin{theorem} \label{EK_n-K_n} 
Under Assumptions \ref{assumptions on kernel 1}, we further assume $\|P\|_\infty = P_{Max}$ for $P_{Max}>0$. Suppose $0<\gamma<\min\{\rho^{-\alpha},1\}$ and $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$. When $\epsilon_n \leq \mathcal{D}_3 \gamma^{\frac{1}{\alpha-d}}$, we have 
\begin{align}\label{Bound variance most general isotropic one}
P\{\|\mathbb{E}K_n-K_n\|_\infty \leq \mathcal{D}_1 \gamma^{1-\frac{d}{\alpha}}\} \geq 1-8(2n)^{2p^2+2p}\exp\left\{-\mathcal{D}_2 n \epsilon_n^{d+\frac{d^2}{\alpha}}\frac{\gamma^{2-\frac{d}{\alpha}}}{N(\gamma)^2}\right\}\,,
\end{align}
where $\mathcal{D}_1$ and $\mathcal{D}_2$ depend on $p$, $d$, $\alpha$, $P_{Max}$, $K_{sup}$, and the second fundamental form of $\iota(M)$, and $\mathcal{D}_3$ depends on $p$, $d$, $\alpha$, $P_{Max}$ and the second fundamental form of $\iota(M)$. 
\end{theorem}

The definitions of $\mathcal{D}_1$, $\mathcal{D}_2$ and $\mathcal{D}_3$ can be found in the proof, which is postponed to Appendix \ref{Appendix A}. The proof involves the methods in \cite{devroye1980strong} developed for the convergence analysis in the Euclidean space setup{, while we need to handle the interaction between the intrinsic and extrinsic geometry}. {We comment that the dependence on the ambient space dimension, that is, $(2n)^{2p^2+2p}$ in the right hand side of \eqref{Bound variance most general isotropic one}, comes from using rectangles in the ambient space to construct a collection of regular sets on the manifold. }


\begin{remark} To appreciate the challenge of generalizing this result to Lebesgue integrable kernels, for simplicity, we assume that $M$ is a rectangle $A \subset \mathbb{R}^d$. Since the function $K(\|\cdot\|_{\mathbb{R}^d})$ is Riemann integrable on $A$, we can approximate $K(\|\cdot\|_{\mathbb{R}^d})$ by a step function $K^*(\cdot)$ uniformly, except on a bad set $A_{\epsilon}$ which can be covered by {\em finitely} many rectangles with the total volume less than $\epsilon>0$. Hence, instead of controlling $|\mathbb{E}K_n(x)-K_n(x)|$ for an arbitrary kernel function, we study the variance of the $0-1$ kernel over rectangles. In contrast, if the function $K(\|\cdot\|_{\mathbb{R}^d})$ is Lebesgue measurable on $A$, we can still approximate $K(\|\cdot\|_{\mathbb{R}^d})$ by a step function $K^*(\cdot)$ uniformly, except on a bad set $B_{\epsilon}$. Here, $B_{\epsilon}$ can be covered by {\em countably}, but may not finitely, many rectangles with the total volume less than $\epsilon$. Therefore, to control the variance $|\mathbb{E}K_n(x)-K_n(x)|$, we have to deal with a bad set $B_{\epsilon}$ associated with each $x$. The major difficulty in generalizing the argument from Riemann integrable kernels to Lebesgue integrable kernels is that if $A_{\epsilon}$ cannot be covered by finite rectangles, then the variance over $B_{\epsilon}$ can not be controlled uniformly for all $x$.  The reader may refer to the last two steps in \eqref{Riemann approx step 3} in the Appendix for details.
\end{remark}


In the following corollary, we state a special case of Theorem \ref{EK_n-K_n} when the kernel $K(t)$ has a compact support. We see that the convergence rate is improved from the general case.

\begin{corollary}\label{compact support variance}
Under Assumption \ref{assumptions on kernel 1}, we further assume that $K(t)$ has a compact support on $[0, \rho]$ and $\|P\|_\infty = P_{Max}$ for $P_{Max}>0$. Suppose $0<\gamma<1$ and $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$. When $\epsilon_n \leq \mathcal{D}_3$, we have
\begin{align}
P\left\{\sup_{x \in M} |\mathbb{E}K_n(x)-K_n(x)| \leq \mathcal{D}_1 \gamma \right\} \geq 1-8(2n)^{2p^2+2p}\exp\left\{-\mathcal{D}_2n \epsilon_n^{d}\frac{ \gamma^2}{N(\gamma)^2}\right\}\,,
\end{align}
where $\mathcal{D}_1$ and $\mathcal{D}_2$ depend on $p$, $d$, $P_{Max}$, $K_{sup}$, and the second fundamental form of $\iota(M)$, and $\mathcal{D}_3$ depends on $\rho$, $P_{Max}$ and the second fundamental form of $\iota(M)$. 
\end{corollary}

The proof of the corollary follows directly by choosing $\gamma<1$ and taking $\alpha \rightarrow \infty$ in Theorem \ref{EK_n-K_n}.  
We explain why the convergence rate is slower when $K(t)$ does not have compact support. As we show in Lemma \ref{rectangle covering lemma}, we have a good control of the variance over any region centered at $x$ that is not too large on the manifold $M$, and this control is uniform. To control the variance outside the region, we have to make sure that the {\em tail} is small enough; that is, for $\eta>0$, $\int_{M \setminus B_\eta(x)} \frac{1}{\epsilon_n^d}K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon_n}\right) P(y)dV(y) $ is small enough when $\epsilon_n$ is small, where $B_\eta(x)$ is a geodesic ball of radius $\eta$ at $x \in M$. When $M$ is an Euclidean space, due to the polynomial decay assumption, the larger $\|y-x\|$ is, the smaller the kernel value $K(\|y-x\|)$ is. So, for $\eta>0$, we can choose $\epsilon_n>0$ sufficiently small so that $\int_{\mathbb{R}^d \setminus B^{\mathbb{R}^d}_{\eta}(x)} \frac{1}{\epsilon_n^d}|K(\frac{\|y-x\|_{\mathbb{R}^d}}{\epsilon_n})| P(y) dy$ is small. 
However, in general, due to the geometry of the manifold, it is not guaranteed that the larger the geodesic distance between $y$ and $x$ is, the smaller the kernel value $K(\|\iota(y)-\iota(x)\|)$. Indeed, two far away points in the sense of geodesic distance might have short Euclidean distance in the ambient space.
Thus, the only way to make sure that $\int_{M \setminus B_\eta(x)} \frac{1}{\epsilon_n^d}K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon_n}\right) P(u) dV(y) $ is small is to choose $\epsilon_n$ {\em further smaller} so that $ \frac{1}{\epsilon_n^d}K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon_n}\right) $ is small outside $B_\eta(x)$. In other words, the slower convergence rate in the non compact kernel case is a result of the geometry of the manifold.

The following Corollary is a special case when the kernel $K(t)$ is a step function. The proof is in  Appendix \ref{Appendix A}.

\begin{corollary}\label{0-1 kernel}
Suppose $K(t)=\sum_{j=1}^J c_j \bigchi_{[a_j, b_j]}(t)$, where $c_j \leq K_{\sup}$ and $0 \leq a_j \leq b_j \leq \rho$. Suppose $\|P\|_\infty = P_{Max}$ for $P_{Max}>0$. Let $0<\gamma<1$ and $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$.  If $\epsilon_n \leq \mathcal{D}_3$, then
\begin{align}
P\left\{\sup_{x \in M} |\mathbb{E}K_n(x)-K_n(x)| \leq \mathcal{D}_1  \gamma\right\} \geq 1-8(2n)^{p+2}\exp\left\{-\mathcal{D}_2  n \epsilon_n^d \gamma^2\right\}\,,
\end{align}
where $\mathcal{D}_1$ depends on $d$, $P_{Max}$, $J$, $K_{\sup}$ and the second fundamental form of $\iota(M)$,  and $\mathcal{D}_2$ depends on $d$, $\rho$, $P_{Max}$ and the second fundamental form of $\iota(M)$ and $\mathcal{D}_3$ depends on $\rho$, $P_{Max}$ and the second fundamental form of $\iota(M)$. 

Hence, if $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$ and $\epsilon_n \leq \mathcal{D}_3$, then with probability greater then $1-\frac{1}{n^2}$, we have 
\begin{align}
\sup_{x \in M} |\mathbb{E}K_n(x)-K_n(x)| \leq \mathcal{D}_4 \sqrt{\frac{\log n}{n \epsilon^d}},
\end{align}
where  $\mathcal{D}_4$ depends on $p$, $d$, $P_{Max}$, $\rho$, $J$, $K_{\sup}$ and the second fundamental form of $\iota(M)$
\end{corollary}

Since $K(t)$ is a step function, we do not need to approximate it anymore. Hence, in contrast to Corollary \ref{compact support variance}, there is no partition number and step function approximation involved and the convergence rate can be further improved. 







\subsection{\textnormal{\textbf{Bias analysis}}}

To get the bias analysis with a convergence rate, we need to further assume the H\"older continuity condition on the probability density function $P$. In fact, if we only need to show that $\|\mathbb{E}K_n(x)-P(x)\|_{L^\infty(M)} \rightarrow 0$ as $\epsilon_n \rightarrow 0$, then it is sufficient to assume $P$ is continuous. However, if we want to control the convergence rate of  $\|\mathbb{E}K_n(x)-P(x)\|_{L^\infty(M)}$, we need to know how $P(x) \rightarrow P(y)$ when $x \rightarrow y$ for $x, y \in M$. 

\begin{Assumption}\label{holder assumption of P}
Suppose the density function $P$ satisfies $\|P\|_\infty = P_{Max}$ for some $P_{Max}>0$ on $M$.  Moreover, $P$ is H\"older continuous so that 
\begin{align}
|P(x)-P(y)| \leq C_P d(x,y)^\kappa\,,
\end{align} 
where $d(x,y)$ is $x,y$ is the geodesic distance between $x$ and $y$ on $M$, $0<\kappa \leq 1$ and $C_P>0$. 
\end{Assumption}

The main theorem in the bias analysis is discussed when $K(t)$ is compactly supported or not. The proof of the theorem is in Appendix \ref{Appendix B}.

\begin{theorem}\label{bias analysis main theorem}
\begin{enumerate}
\item Under Assumptions \ref{assumptions on kernel 1} and \ref{holder assumption of P}, when $0<\gamma<1$ and $\epsilon \leq \omega_1 \gamma^{\frac{2\alpha}{(\alpha-d)}}$, we have
\begin{align}
\sup_{x \in M}\left|\int_{M}\frac{1}{\epsilon^d} K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon}\right)P(y)dV(y)-P(x)\right|\leq \omega_2 \gamma^\kappa\,,
\end{align}
where $\omega_1$ depends $\rho$, $\alpha$ and $d$ and $\omega_2$ depends on $\rho$, $\alpha$, $d$, $\kappa$, $C_P$, $P_{Max}$, $K_{sup}$, the curvature of $M$ and the second fundamental form of $\iota(M)$.

\item Under Assumptions \ref{assumptions on kernel 1} and \ref{holder assumption of P}, we further assume that $K(t)$ is compactly support on $[0, \rho]$, where $\rho>0$. When $0<\gamma<1$ and $\epsilon \leq \omega_1 \gamma^{2}$, we have
\begin{align}
\sup_{x \in M}\left|\int_{M}\frac{1}{\epsilon^d} K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon}\right)P(y)dV(y)-P(x)\right|\leq \omega_2 \gamma^\kappa\,,
\end{align}
where $\omega_1$   depends on $\rho$ and $\omega_2$  depends on $\rho$, $d$, $\kappa$, $C_P$, $P_{Max}$, $K_{sup}$, the curvature of $M$ and the second fundamental form of $\iota(M)$.
\end{enumerate}
\end{theorem}




\subsection{\textnormal{\textbf{Put everything together}}}

We concatenate the variance and the bias analysis to derive the main result. 
In the following theorem, we discuss three cases corresponding to those three cases in the variance analysis. 
The proof is in Appendix \ref{Appendix B}.

\begin{theorem}\label{density estimation}
\begin{enumerate}
\item Under Assumptions \ref{assumptions on kernel 1} and \ref{holder assumption of P}, assume $0<\gamma<\min\{\rho^{-\alpha},1\}$ and $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$. When $\epsilon_n \leq \Omega_3 \gamma^{\frac{2\alpha}{\alpha-d}}$, we have 
\begin{align}
&P\left\{\sup_{x \in M}|K_n(x)-P(x)| \leq \Omega_1 (\gamma^{1-\frac{d}{\alpha}}+\gamma^\kappa)\right\} \nonumber \\
\geq \,&1-8(2n)^{2p^2+2p}\exp\left\{-\Omega_2 n \epsilon_n^{d+\frac{d^2}{\alpha}}\left(\frac{\gamma^{2-\frac{d}{\alpha}}}{N(\gamma)^2}\right)\right\}\,,
\end{align}
where $\Omega_1$ depends on $p$, $\rho$, $d$, $\alpha$, $\kappa$, $C_P$, $P_{Max}$, $K_{sup}$, the curvature of $M$ and the second fundamental form of $\iota(M)$. $\Omega_2$  depends on $p$, $d$, $\alpha$, $P_{Max}$, $K_{sup}$ and the second fundamental form of $\iota(M)$. $\Omega_3$  depends on $p$, $\rho$, $\alpha$, $d$, $P_{Max}$ and the second fundamental form of $\iota(M)$. 

\item Under Assumptions \ref{assumptions on kernel 1} and \ref{holder assumption of P}, we further assume that $K(t)$ is compactly supported on $[0, \rho]$. Suppose $0<\gamma<1$ and $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$.  If $\epsilon_n \leq \Omega_3\gamma^{2}$, then
\begin{align}
P\left\{\sup_{x \in M}|K_n(x)-P(x)| \leq \Omega_1 \gamma^\kappa \right\} \geq 1-8(2n)^{2p^2+2p}\exp\left\{-\Omega_2\frac{n \epsilon_n^{d} \gamma^2}{N(\gamma)^2}\right\}\,,
\end{align}
where $\Omega_1$  depends on $p$, $\rho$, $d$, $\kappa$, $C_P$, $P_{Max}$, $K_{sup}$, the curvature of $M$ and the second fundamental form of $\iota(M)$. $\Omega_2$  depends on $p$, $d$, $P_{Max}$, $K_{sup}$, and the second fundamental form of $\iota(M)$. $\Omega_3$  depends on $\rho$,  $P_{Max}$ and the second fundamental form of $\iota(M)$. 
\item
Suppose $K(t)=\sum_{j=1}^J c_j \bigchi_{[a_j, b_j]}(t)$, where $c_j \leq K_{\sup}$ and $0 \leq a_j \leq b_j \leq \rho$. Moreover, $K(t)$ satiesfies (iv) in Assumption \ref{assumptions on kernel 1} . Suppose $0<\gamma<1$ and $\epsilon_n \rightarrow 0$, as $n \rightarrow \infty$. Under Assumption \ref{holder assumption of P}, if $\epsilon_n \leq \Omega_3\gamma^{2}$, then
\begin{align}
P\left\{\sup_{x \in M}\left|K_n(x)-P(x)\right| \leq \Omega_1 \gamma^\kappa \right\} \geq 1-8(2n)^{p+2}\exp\left\{-\Omega_2  n \epsilon_n^d \gamma^2\right\}\,,
\end{align}
where $\Omega_1$ depends on $\rho$, $d$, $J$, $\kappa$, $C_P$, $P_{Max}$, $K_{\sup}$, the curvature of $M$ and the second fundamental form of $\iota(M)$. $\Omega_2$ which depends on $d$, $\rho$, $P_{Max}$ and the second fundamental form of $\iota(M)$. $\Omega_3$ which depends on $\rho$,  $P_{Max}$ and the second fundamental form of $\iota(M)$. 

Hence, if $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$ and $\epsilon_n \leq \Omega_3$, then with probability greater then $1-\frac{1}{n^2}$, we have 
\begin{align}
\sup_{x \in M} |K_n(x)-P(x)| \leq \Omega_4 (\frac{\log n}{n \epsilon^d})^{\frac{\kappa}{2}},
\end{align}
where  $\Omega_4$ depends on $p$, $\rho$, $d$, $J$, $\kappa$, $C_P$, $P_{Max}$, $K_{\sup}$, the curvature of $M$ and the second fundamental form of $\iota(M)$.
\end{enumerate}
\end{theorem}

\subsection{Riemann integrability of an isotropic kernel on a manifold}\label{Relationship between intrinsic and extrinsic}

In this section, {we discuss the delicate difference when the kernel is defined on the ambient Euclidean space and on the manifold, particularly when the kernel is irregular. The interaction between intrinsic and extrinsic geometry plays a role.}
We show a necessary and sufficient condition for an isotropic kernel to be Riemann integrable on a manifold when it is Riemann integrable in the ambient space. We will provide an {explicit example of a kernel which} is Riemann integrable in the ambient space but not integrable on the manifold. {In the example}, although the kernel is not Riemann integrable on the manifold, by our main theorem, such kernel can still be used for density estimation on the manifold and the convergence rate can still be estimated. Our theorem in this section {suggests} that the results in the Euclidean space may not be easily generalized to the manifold case when the extrinsic properties of the manifold are involved. 


For any function $f$ which is Riemann integrable on any closed ball in $\mathbb{R}^p$, it is not necessary that $f$ is a Riemann integrable function over any embedded submanifold in $\mathbb{R}^p$. A trivial example is as follows.

\begin{example}
$f(x,y)$ is defined on $\mathbb{R}^2$. $f(x,y)=1$ when $(x,y) \in \mathbb{Q}\cap (0,1) \times \{0\}$. And $f(x,y)=0$ otherwise. Then $f(x,y)$ is Riemann integrable over any  ball in $\mathbb{R}^p$. But, if $(0,1) \times \{0\}$ is the embedded submanifold, then $f(x,y)$ is not Riemann integrable over the submanifold.
\end{example}

For the theoretical purpose, one may ask that whether $K(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon})$ is a Riemann integrable function on the manifold $M$. 
The next theorem provides the necessary and sufficient condition such that $K(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon})$ is Riemann integrable on the manifold $M$ for any $\epsilon$. 

\begin{theorem}\label{Riemann integrable main theorem}
Suppose $M$ is a $d$-dimensional compact smooth manifold without boundary isometrically embedded in $\mathbb{R}^p$ through $\iota$. Fix $x \in M$. Let $D_x(\textbf{u})=\|\textbf{u}-\iota(x)\|_{\mathbb{R}^p}$. The set of critical points of $D_x(\textbf{u})$ on $\iota(M)$ is Jordan measurable if and only if for any $\epsilon>0$, $K(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon})$ is a Riemann integrable function of $y$ on the manifold for all bounded kernel $K(t):\mathbb{R}_{\geq0}\to \mathbb{R}$ that is Riemann integrable on $[0,a]$ for all $a>0$.
\end{theorem}

The proof of Theorem \ref{Riemann integrable main theorem} is in Appendix \ref{Appendix D}. When $M$ is an analytic manifold, the set of critical points of $D_x(\textbf{u})$ on $\iota(M)$ has measure $0$. Because the set of critical points of $D_x(\textbf{u})$ on $\iota(M)$ is a closed subset of $\iota(M)$ containing $\iota(x)$, it contains all its boundary points. Hence, the set of critical points of $D_x(\textbf{u})$ on $\iota(M)$ is Jordan measurable and we have the following corollary.

\begin{corollary}
Suppose $M$ is a $d$-dimensional compact {\em analytic} manifold without boundary isometrically embedded in $\mathbb{R}^p$ through $\iota$. Take a bounded kernel $K(t):\mathbb{R}_{\geq0}\to \mathbb{R}$ that is Riemann integrable on $[0,a]$ for all $a>0$.  Fix $x \in M$. For any $\epsilon>0$, $K(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon})$ is a Riemann integrable function of $y$ on the manifold.
\end{corollary}

Next, we construct explicitly an embedded manifold $\iota(M)$ so that for some $\iota(x)$ on $\iota(M)$, the set of critical points of $D_x(\textbf{u})$ on $\iota(M)$ is not Jordan measurable. Moreover, we construct a function $K(t)$ satisfies Assumption \ref{assumptions on kernel 1}, but $K(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon})$  is not Riemann integrable on the manifold for infinity many choices of $\epsilon$.

\begin{example} 
We construct a fat Cantor set $C$ in $[0, \frac{\pi}{2}]$ with non-zero measure. Then there is a non-negative smooth function $f(\theta)$ that only vanishes on $C$.  $f(\theta)$ can be constructed in the following way. The complement of $C$ in $[0, \frac{\pi}{2}]$ is the union of countable open intervals. Then on each of the open interval, $f(\theta)$ is equal to an everywhere positive bump function subject to  the open interval and vanishes on the boundary of the open interval. Moreover, $f(\theta)$ is $0$ on the fat cantor set. We construct  the curve $\gamma(\theta) \subset \mathbb{R}^2$ for $\theta \in  [0, 2\pi]$ with the following conditions. 
\begin{enumerate}
\item $\gamma(\theta)$ is smooth with $\gamma(\pi)=(0,0)$.
\item $\gamma(\theta)=((f(\theta)+1)cos(\theta), (f(\theta)+1)sin(\theta))$, for $\theta \in [0, \frac{\pi}{2}]$.  
\item $\gamma(\theta)$ for $\theta \in (\frac{\pi}{2}, 2\pi)$ is contained in the open unit disc centered at $(0,0)$.
\end{enumerate}
Note that the fat Cantor set $C$ is embedded into the closed curve as a fat Cantor set. In other words, let $C'=\gamma(C)$. Then, the unit circle in $\mathbb{R}^2$ is tangent to $\gamma(\theta)$ along the fat cantor set $C'$.  Hence, they are the critical points of the function $D_x(\textbf{u})=\|\textbf{u}\|_{\mathbb{R}^2}$. The other critical points on $\gamma(\theta)$ for $\theta \in [0, \frac{\pi}{2}]$ are the image of the critical points of $f(\theta)$ in the open intervals under $\gamma$. Thus, they are countable. Therefore, the set of critical points of $D_x(\textbf{u})=\|\textbf{u}\|_{\mathbb{R}^2}$ is not Jordan measurable.
\begin{align}
K(t)=\left\{
\begin{array}{lll} 
\frac{1}{3}& \mbox{if $t \in [0, 1)\cup (1,\frac{3}{2}]$};\\ \frac{1}{t^2} & \mbox{if $t$ is a positive integer};\\ 
0 & \mbox{everywhere else}. \end{array} \right.
\end{align}
Note that $K(t)$ satisfies Assumption \ref{assumptions on kernel 1} and $K(t)$ is not continuous on positive integers and $\frac{3}{2}$. Fix $x=0\in \mathbb{R}^2$. We can easily see that $K(\frac{\|y-0\|_{\mathbb{R}^2}}{\epsilon})=K(\frac{\|y\|_{\mathbb{R}^2}}{\epsilon})$ is not a Riemann integrable function  on $\gamma(\theta)$, for any $\epsilon=\frac{1}{k}$ where $k$ can be any positive integer.   In other word, given the above manifold and kernel, there are arbitrarily small bad choices of $\epsilon$ so that the kernel fails to be Riemann integrable. 
\end{example}

 \subsection{A comparison between Riemann integrable kernels and kernels in the VC class}

We recall the VC class of kernels subject to a set $A$ in $\mathbb{R}^p$, and refer readers with interest to \cite{shorack2009empirical} and \cite{pollard2012convergence} for more general definition. Suppose $K(\textbf{u})$ is a bounded function on $\mathbb{R}^p$ and $K(\textbf{u}) \in L^1 (\mathbb{R}^p)$. Note that $K(\textbf{u}) $ is also in  $ L^2 (\mathbb{R}^p)$.   For any $A\subset\mathbb{R}^p$, we consider the {\em space of kernels} over $A$, denoted as $\mathcal{F}(A)$, by 
\begin{align}
\mathcal{F}(A)=\{K(x - \cdot ), \, x\in A \}.
\end{align}
Suppose $\mathcal{P}$ is any probability measure defined on the $\sigma$-algebra of the Borel sets in $\mathbb{R}^p$. Then the $L^2(\mathcal{P})$ metric over $\mathcal{F}(A)$ is defined as 
\begin{align}\label{general L2 metric}
d_{L^2(\mathcal{P})}(K(x - \cdot ), K(y - \cdot ))=(\int_{\mathbb{R}^p} (K(x - z)-K(y -z))^2 d \mathcal{P}(z) ) ^{\frac{1}{2}}.
\end{align}
Let $N_{cov}(\epsilon, \mathcal{F}(A), d_{L^2(\mathcal{P})})$ be the $\epsilon$-covering number of $\mathcal{F}(A)$ with respect to the metric $d_{L^2(\mathcal{P})}$. Then $\mathcal{F}(A)$ is a VC class, whenever there exist constants $C>0$ and $b>0$ such that for all $0<\epsilon<1$,
\begin{align}
\sup_{\mathcal{P}} N_{cov}(\epsilon, \mathcal{F}(A), d_{L^2(\mathcal{P})}) \leq C \epsilon^{-b},
\end{align}
where the supremum is taken over all the probability measure $\mathcal{P}$ defined on the $\sigma$ algebra of the Borel sets in $\mathbb{R}^p$.
The constants $C$ and $b$ are called the VC characteristics. In \cite{shorack2009empirical}, \cite{pollard2012convergence}, \cite{nolan1987u}, \cite{van2000applications}, the authors discuss several sufficient conditions for $\mathcal{F}(A)$ to be a VC class. For example, if $g_1$ is a bounded real function with bounded  variation, $g_2(\textbf{u})$ is a polynomial on $\mathbb{R}^p$ and $K(\textbf{u})=g_1(g_2(\textbf{u}))$, then $\mathcal{F}(A)$ is a VC class.

Next, we show an example of kernel which is Riemann integrable but not in VC class. 
\begin{example} \label{not VC class}
Let $K(t)=\sin(\exp(\exp (\frac{1}{|t|})))$ for $|t| \leq 1$ and $K(t)=0$ for $|t|>1$. $K(t)$ is discontinuous at $t=-1, 0, 1$. Hence, $K(t)$ is Riemann integrable. It is also trivial to modify $K(t)$ so it is discontinuous only at  $t=0$. Suppose $A=[0,1] \subset \mathbb{R}$, then $\mathcal{F}(A)$ is not a VC class. The proof is in Appendix \ref{proof of not VC class}. 
\end{example}

\begin{example} \label{not RI class}
There are kernels in the VC class that are not Riemann integrable. For example, suppose $K_1$ is the characteristic function on the irrational numbers in $[0,1]$. $K_1$ is clearly in the VC class but not Riemann integrable. 
To further explore $K_1$, suppose $K_2$ is the characteristic function on $[0,1]$. Clearly, $K_2$ differs from $K_1$ by a measure $0$ set, and $K_1$ and $K_2$ have the same $\epsilon$ covering number of any closed interval $A$ on $\mathbb{R}$.  However,  it can be shown that $K_1$ and $K_2$ have the same $L^\infty$ convergence rate when they are applied for KDE. In other words, the KDE behavior of the non-Riemann integrable kernel $K_1$ can be studied through the Riemann integrable kernel $K_2$. 
\end{example}


Before closing this subsection, we have a comparison of our results with those shown in \cite{gine2002rates}. The authors in \cite{gine2002rates} proved the following variance analysis result for the density estimation on the Euclidean space $\mathbb{R}^d$ by using the kernel in the VC class.
\begin{theorem}\label{Gine Guillou}[Gin\'e and Guillou]
Let $P$ be a bounded probability density function on $\mathbb{R}^d$. Let $\textbf{v}_1, \cdots, \textbf{v}_n$ be a sequence of i.i.d samples from $\mathbb{R}^d$ based on $P$. Suppose $K(\textbf{v})$ is in the VC class of $\mathbb{R}^d$ with the VC characteristics $C$ and $b$. Suppose
\begin{align}\label{condition 1 in Gine}
 \|P\|_{\infty}  \int_{\mathbb{R}^d} |K(\textbf{v})|^2d\textbf{v} \leq D.
\end{align}
For $\textbf{v} \in \mathbb{R}^d$, define the kernel density estimator at $\textbf{v}$ to be
\begin{align}
K_{n,\epsilon}(\textbf{v})=\frac{1}{n \epsilon^d}\sum_{i=1}^n K\left(\frac{\textbf{v}_i-\textbf{v}}{\epsilon}\right)\,. 
\end{align}
Suppose $C_1$ and $C_2$ are constants depending on the VC characteristics. For any $c_1>C_1$ and $0<\gamma<\frac{c_1D}{\|K\|_\infty}$, there is a $n_0$ depending on $\gamma$, $D$, $\|K\|_\infty$ and the VC characteristics, such that if $n>n_0$, then
\begin{align}\label{GiniMainresult}
P\left\{\sup_{x \in \mathbb{R}^d} |\mathbb{E}K_n(\textbf{v})-K_n(\textbf{v})| \geq 2 \gamma\right\} \leq C_2 \exp\left\{-\frac{1}{D}\frac{\log(1+\frac{4c_1}{C_2})}{c_1C_2}  n \epsilon^d \gamma^2\right\}\,,
\end{align}
\end{theorem}

It is reasonable to compare Theorem \ref{Gine Guillou} with  Corollary \ref{compact support variance}. In Corollary \ref{compact support variance}, we prove that if $\epsilon_n \leq \mathcal{D}_3$, then
\begin{align}\label{OurCorollary1Statement}
P\left\{\sup_{x \in M} |\mathbb{E}K_n(x)-K_n(x)| \geq \mathcal{D}_1 \gamma \right\} \leq 8(2n)^{2p}\exp\left\{-\mathcal{D}_2n \epsilon_n^{d}\frac{ \gamma^2}{N(\gamma)^2}\right\}\,,
\end{align}
which has an extra $\log(n)$ term compared with \eqref{GiniMainresult}.
The term $\frac{\log(1+\frac{4c_1}{C_2})}{c_1C_2}$ in \eqref{GiniMainresult} and the term $\frac{1}{N(\gamma)^2}$ in \eqref{OurCorollary1Statement} both characterize the regularity of the kernel. Besides this difference, the convergence rate in Theorem \ref{Gine Guillou} and Corollary \ref{compact support variance} are the same.  However, in the case when the kernel is not in the VC class,  $\frac{\log(1+\frac{4c_1}{C_2})}{c_1C_2}$ does not exist.  
{
It would be an interesting future direction to generalize the $L^\infty$ convergence result with a kernel in the VC class and provide the rate in the manifold setup. 

}


{
\section{Kernel density estimation with non-translational invariant Riemann integrable kernels}\label{Section:KDE3} 

We switch our discussion to more general Riemann integrable kernels that might be non-translational invariant. As we discussed in the introduction, this is the case that we may encounter when we analyze modern machine learning algorithms. Note that when the kernel function is defined on the Euclidean space, the kernel's property might be different when it is restricted on the manifold, as is shown in Section \ref{Relationship between intrinsic and extrinsic}. In this case, we thus need to further know the kernel's behavior on the manifold so that we can guarantee if the density estimator works.
In the field of manifold learning, sometimes the kernel is implicitly defined on the manifold via the algorithm. In this case, we could focus on the kernel analysis. 

We focus on a special non-translational invariant kernel that is inspired by the LLE algorithm. 
Suppose a manifold learning algorithm implicitly provides us with kernel via a chart of $\iota(M)$, where the kernel might be irregular. We can use this kernel to construct a kernel density estimator. 
We need the following two assumptions to describe the general kernel function.  

\begin{Assumption}\label{assumptions on non isotropic kernel}
\begin{enumerate}[(i)]
\item
$K(\textbf{v}):\mathbb{R}^d \rightarrow \mathbb{R}$ is a bounded function such that $\sup|K(\textbf{v})| = K_{\sup}$ for some $K_{\sup}>0$.
\item
$K$ is Riemann integrable on $\mathbb{R}^d$ with the support contained in the cube $[-R,R]^d$.
\item
$\int_{ \mathbb{R}^d} K(\textbf{v})d\textbf{v}=1$.
\end{enumerate}
\end{Assumption}

Next, suppose the chart satisfies the following assumption.  Recall that with a chart, for each point $x\in M$, there is a local diffeormorphism from $M$ to some Euclidean space.
\begin{Assumption}\label{assumptions on charts of M}
Suppose there are constants $r$, $D_1 \geq 1$, and $D_2\geq0$ such that for any $x \in M$ there is a diffeomorphism $\Phi_x$ satisfying the following conditions: 
\begin{enumerate}[(i)]
\item
$\Phi_x(\textbf{v}): B^{\mathbb{R}^d}_r(0) \rightarrow \iota(M)$ with $\Phi_x(0)=\iota(x)$.  $B^{\mathbb{R}^d}_r(0)$ is the open Euclidean ball of radius $r$ centered at $0$.
\item
For $\textbf{v}_1, \textbf{v}_2 \in B^{\mathbb{R}^d}_r(0)$, we have  $\frac{1}{D_1} \|\textbf{v}_1- \textbf{v}_2\|_{\mathbb{R}^d} \leq d_{\iota(M)} (\Phi_x(\textbf{v}_1),\Phi_x(\textbf{v}_2))  \leq  D_1 \|\textbf{v}_1- \textbf{v}_2\|_{\mathbb{R}^d} $, where $d_{\iota(M)}$ is the geodesic distance on $\iota(M)$. 
\item
Suppose the volume form of $\iota(M)$ in the chart $\Phi_x$ can be expressed as $dV=U_x(\textbf{v}) d \textbf{v}$ for $\textbf{v} \in B^{\mathbb{R}^d}_r(0)$.  Then, $|\frac{U_x(\textbf{v})}{U_x(0)}-1| \leq D_2 \|\textbf{v}\|_{\mathbb{R}^d}$.
\end{enumerate}
\end{Assumption}
Note that (ii) in the above assumption implies that $|U_x(\textbf{v})|\leq D^d_1$ for any $\textbf{v} \in B^{\mathbb{R}^d}_r(0)$. Hence, (iii) is equivalent to the fact that we have a $C^1$ bound of the volume form at $0$. 
Now, we define the kernel density estimator. For $x \in M$, define the kernel density estimator at $x$ with the bandwidth $\epsilon_n>0$ as
\begin{align}
K_{n}(x)=\frac{1}{n \epsilon_n^d U_x(0)}\sum_{i=1}^n K_{\epsilon_n}(x,x_i)\,,
\end{align}
where 
\[
K_{\epsilon_n}(x,x_i):=K\left(\frac{\Phi^{-1}_x(\iota(x_i))}{\epsilon_n}\right)
\] 
if $\iota(x_i) \in \Phi_x(B^{\mathbb{R}^d}_r(0))$ and $K_{\epsilon_n}(x,x_i)=0$ otherwise. Note that the kernel function $K$ is defined on $\mathbb{R}^d$, which can be understood as the tangent space. Also note that the kernel is not translational invariant since in general $K_{\epsilon_n}(x,x_i)$ cannot be written as $k(x-x_i)$ for some function $k$ defined on $\mathbb{R}^p$. In other words, at different points, we use different weights to estimate the density function. Clearly, this is a special non-translational invariant kernel. In general, even if we have the chart, we still need to estimate $U_x(0)$. In the special case when $\Phi_x (\textbf{v})= \iota \circ \exp_x(\textbf{v})$ where $\exp_x$ is the exponential map at $x$, an expansion of $U_x(\textbf{v})$ in terms of $\textbf{v}$ is shown in Lemma \ref{volume form in normal coordinate} in Appendix. Especially,  $U_x(0)=1$ in this case.


Now we claim the strong uniform consistency of this kernel density estimator. Define the partition number $N(\gamma)$ of $K(\textbf{v})$ over $[-R, R]^d$ in the same way as the partition number for an isotropic kernel in \eqref{equation def partition number}. We have the following result.
\begin{theorem}\label{Theorem nontranslational invariant Riemann}
Suppose Assumptions  \ref{holder assumption of P}, \ref{assumptions on non isotropic kernel}, and \ref{assumptions on charts of M} hold. Take $0<\gamma<1$ and assume $\epsilon_n \rightarrow 0$ as $n \rightarrow \infty$.  If $\epsilon_n \leq \Omega_3\gamma^{2}$, then
\begin{align}
P\left\{\sup_{x \in M}|K_n(x)-P(x)| \leq \Omega_1 \gamma^\kappa \right\} \geq 1-8(2n)^{2p^2+2p}\exp\left\{-\Omega_2\frac{n \epsilon_n^{d} \gamma^2}{N(\gamma)^2}\right\}\,,
\end{align}
where $\Omega_1$  depends on $p$, $R$, $d$, $\kappa$, $C_P$, $P_{Max}$, $K_{sup}$, $D_1$, $D_2$, and the second fundamental form of $\iota(M)$. $\Omega_2$  depends on $p$, $d$, $P_{Max}$, $K_{sup}$, $D_1$, and the second fundamental form of $\iota(M)$. $\Omega_3$  depends on $r$, $R$, $D_1$, $P_{Max}$, and the second fundamental form of $\iota(M)$. 
\end{theorem}

Thanks to the chart, the proof of the above theorem is a direct extension of the variance and the bias analysis for the isotropic kernel case shown in Section \ref{Section:KDE2}. Since the proof is similar, we omit details but only indicate key steps. Specifically, by (ii)  in Assumption \ref{assumptions on charts of M} and the argument in Lemma \ref{approximation to K} in Appendix, we can find an approximation of $K(\frac{\Phi^{-1}_x(\iota(y))}{\epsilon_n})$ by a simple function. The rest of the proof in the variance analysis is same as the proof of Theorem \ref{EK_n-K_n}. The bias analysis follows the same argument as in the proof of Theorem \ref{bias analysis main theorem} by applying the change of variable formula in the chart $\Phi_x$. 



The kernel function associated with LLE \cite{wu2018think} that motivates the above result deserves a discussion. 
Take $\rho\in\mathbb{R}$ as the regularization order in the LLE algorithm \cite[See (2.9) and (3.5)]{wu2018think} and $\epsilon>0$ as a bandwidth. The kernel associated with {LLE} is  
\begin{equation}
K_{\texttt{LLE}}(x,y)=[1- \mathbf{T}_{\iota(x)}^\top(\iota(y)-\iota(x))]\chi_{B_{\epsilon}^{\mathbb{R}^p}(\iota(x)) \cap \iota(M)}(\iota(y)),
\end{equation}
where $x,y\in M$, 
\begin{equation}
\mathbf{T}_{\iota(x)}:= \mathcal{I}_{\epsilon^{d+\rho}}(C_{x})\big[\mathbb{E}(X-\iota(x))\chi_{B_{\epsilon}^{\mathbb{R}^p}(x)}\big]    \in \mathbb{R}^p\,,\label{Definition:Tx:ContinuousCase}
\end{equation}
\[
C_x:=\mathbb{E}[(X-\iota(x))(X-\iota(x))^{\top}\chi_{B_{\epsilon}^{\mathbb{R}^p}(\iota(x))}(X)]\,
\]
and
$\mathcal{I}_{\epsilon^{d+\rho}}(C_{x})$ is the regularized pseudoinverse of $C_x$ with the regularization $\epsilon^{d+\rho}$.   
Here, $C_x$ is understood as the {\em local covariance matrix}, which depends on the embedding and hence the extrinsic geometry of the manifold. This quantify is commonly considered in the local principal component analysis. Moreover, $\mathbf{T}_{\iota(x)}$ depends on the manifold curvature, mainly the extrinsic curvature information.


Note that the LLE algorithm is based on the barycentric coordinate, and its relationship with KDE by using $K_{\texttt{LLE}}$ is not obvious by reading the algorithm. In \cite[Sections 2 and 3]{wu2018think}, the barycentric coordinate is carefully rewritten so that locally points are represented as coordinates associated with a section of the associated frame bundle. Such conversion naturally leads to the KDE analysis. The density estimation by using $K_{\texttt{LLE}}$ appears naturally in the algorithm as is shown in \cite[(3.15)]{wu2018think}.
Clearly, in general the kernel function $K_{\texttt{LLE}}$ is non-translational invariant, and it is not designed by the user but emerges implicitly in the algorithm. Moreover, when the scaling parameter $\epsilon$ changes and the density function changes, the kernel function itself changes correspondingly. Next, we consider a special case.

\begin{example}
Suppose $M=S^{p-1}$, the canonical $(p-1)$-dim sphere, and it is isometrically embedded in $\mathbb{R}^p$ through $\iota$ with the center at the origin. To simplify the discussion, we also assume the data is uniformly sampled. By a direct calculation in \cite[Appendix G.1]{wu2018think} and the cosine law, we can describe $K_{\texttt{LLE}}$ as follows. For $\textbf{v} \in \mathbb{R}^{p-1} $, let 
\begin{align}
K_{\epsilon}(\textbf{v}) = [1-a(\epsilon, \rho)+a(\epsilon, \rho)\cos(\|\textbf{v}\|_{\mathbb{R}^{p-1}})]\chi_{B_{r(\epsilon)}^{\mathbb{R}^{p-1}}}(\textbf{v}), 
\end{align}
where $r(\epsilon)=\arccos(1-\frac{\epsilon^2}{2})$ and $a>0$ is a function depending on $\epsilon$ and $\rho$.  Then, $K_{\texttt{LLE}}(x,y)=K_{\epsilon}(\exp^{-1}_x(y))$ for $y$ in the geodesic ball of radius $r(\epsilon)$ around $x$ and $K_{\texttt{LLE}}(x,y)=0$ otherwise. 
Based on the above discussion, we see that $K_{\texttt{LLE}}(x,y)$ is a kernel depending on the intrinsic structure (the normal coordinates) of the manifold. Due to the symmetry of $S^{p-1}$, $K_{\texttt{LLE}}(x,y)$ can be further simplified. In this special case, $\|\exp^{-1}_x(y)\|_{\mathbb{R}^{p-1}}=\arccos(1-\frac{\|\iota(x)-\iota(y)\|_{\mathbb{R}^p}^2}{2})$. Hence, if we take
\begin{align}
\tilde{K}_{\epsilon}(t) = \left[1-\frac{a(\epsilon, \rho)\epsilon^2}{2} t^2\right]\chi_{[0,\epsilon]}(t), 
\end{align}
$K_{\texttt{LLE}}(x,y)=\tilde{K}_{\epsilon}(\frac{\|\iota(x)-\iota(y)\|_{\mathbb{R}^p}}{\epsilon})$.
When $\rho \leq 4$, there is a constant $C$ such that $0<a(\epsilon, \rho)\epsilon^2<C$ for all $\epsilon$ small enough. Therefore, if $\rho \leq 4$ and $\epsilon$ is small enough, $K_{\epsilon}(\|\textbf{v}\|_{\mathbb{R}^p})$ has a uniform $C^1$ bound except for $\textbf{v}$ on the boundary of $B_{\epsilon}^{\mathbb{R}^p}(0)$ and $K_{\epsilon}(\|\textbf{v}\|_{\mathbb{R}^p})$ has a $L^\infty$ bound for $\textbf{v}$ on the boundary of $B_{\epsilon}^{\mathbb{R}^p}(0)$. The volume of the boundary of $B_{\epsilon}^{\mathbb{R}^p}(0)$ becomes smaller when $\epsilon$ is smaller. Hence, there is a uniform partition so that \eqref{equation def partition number} holds for all $K_{\epsilon}(\|\textbf{v}\|_{\mathbb{R}^p})$ with $\epsilon$ small enough. We conclude that we have an upper bound of the partition numbers for all  $K_{\epsilon}(\|\textbf{v}\|_{\mathbb{R}^p})$ that only depends on $\gamma$ and is independent of $\epsilon$. Hence, the $L^\infty$ convergence and its rate of the density estimator $\frac{p-1}{|S^{p-2}| n \epsilon_n^{p-1}} \sum_{i=1}^n K_{\texttt{LLE}}(x, x_i)$ for $\rho \leq 4$ are the same as those shown in Part 2, Theorem \ref{density estimation} and can be proved by the same method.
\end{example}

 In this simple case, due to the symmetry, the kernel is isotropic. In general, the LLE kernel is far from being so nice, and its strong uniform consistency cannot be inferred by  Theorem \ref{density estimation} or Theorem \ref{Theorem nontranslational invariant Riemann}. While the LLE kernel is in general not too irregular, the main challenge is that the kernel, and hence the partition number, depends on $\epsilon$.
In general, the convergence analysis should depend on the algorithm, where the kernel might be not only non-translational invariant, but also dependent on $\epsilon$ and irregular. We will report such analyses in our future research.
}




\section{Kernel density estimation with Lebesgue integrable kernels}\label{Section:KDE4} 


We consider the following assumption on the Lebesgue integrable kernel.

\begin{Assumption}\label{assumptions on kernel}
$K_{\epsilon}:\iota(M) \times \iota(M) \rightarrow \mathbb{R}_{\geq 0}$, where $\epsilon>0$, is a sequence of kernel functions so that the following conditions are satisfied for all $\epsilon$:
\begin{enumerate}[(i)]
\item There exists a constant $K_{\sup}>0$, such that $0 \leq K_\epsilon(\iota(x),\iota(y)) \leq \frac{K_{\sup}}{\epsilon^\alpha}$, for all $x,y \in M$ and $\alpha \geq d$.

\item $K_{\epsilon}$ is a Lebesgue measurable function on $\iota(M) \times \iota(M)$.

\item $K_{\epsilon}(\iota(x),\iota(y))=0$ if $\|\iota(x)-\iota(y)\|_{\mathbb{R}^p}>\epsilon$.

\item $\int_{M} K_{\epsilon}(\iota(x),\iota(y))dV(y)=1$, where $dV$ is the volume form on $M$.
\end{enumerate}
\end{Assumption} 

\begin{remark}
For the kernel $K(t)$ satisfies Assumption \ref{assumptions on kernel 1}, if we apply the same notation as in Assumption  \ref{assumptions on kernel}, we have $K_{\epsilon}(\iota(x), \iota(y))=\frac{1}{\epsilon^d}K\left(\frac{\|\iota(y)-\iota(x)\|_{\mathbb{R}^p}}{\epsilon}\right)$. Hence, for the isotropic kernel,  $\epsilon^d K_{\epsilon}(\iota(x), \iota(y)) \leq K_{\sup}$  for all $\epsilon$. In contrast, if $K_{\epsilon}$ is a kernel defined as in  Assumption  \ref{assumptions on kernel}, then $\epsilon^d K_{\epsilon}(\iota(x), \iota(y)) \leq \frac{K_{\sup}}{\epsilon^{\alpha-d}} $. In other words, Assumption  \ref{assumptions on kernel} is more general than Assumption \ref{assumptions on kernel 1} in the sense that we allow the kernel to blow up as $\epsilon \rightarrow 0$.
\end{remark}

Let $x_1, \cdots, x_n$ be a sequence of i.i.d samples from $M$ based on the probability density function $P$ on $M$.  For $x \in M$, we define the kernel density estimator at $x$ to be
\begin{align}
K_{n,\epsilon}(x)=\frac{1}{n}\sum_{i=1}^n K_{\epsilon}(\iota(x),\iota(x_i))\,.
\end{align}
Clearly, we have
\begin{align}
\mathbb{E}K_{\epsilon}(x)=\int_{M}K_{\epsilon}(\iota(x),\iota(y))P(y)dV(y).
\end{align}
Then we have the following result for the variance analysis. 
\begin{proposition}\label{variance L1 case}
Suppose $P$ is measurable on $M$ with $\|P\|_\infty = P_{max}$ for some $P_{max}>0$. Under Assumption \ref{assumptions on kernel}, if $\epsilon \rightarrow 0$ and  $\frac{\log n}{n \epsilon^{2\alpha-d}} \rightarrow 0$ as $n \rightarrow \infty$, then, a.s.
\begin{align}
\int_{M} \left|\frac{1}{n}\sum_{i=1}^n K_{\epsilon}(\iota(x),\iota(x_i))-\mathbb{E}K_{\epsilon}(x)\right| dV(x) \rightarrow 0. 
\end{align}
\end{proposition}

The proof of the proposition is in Appendix \ref{Appendix C}, which generalizes the method developed in \cite{devroye1979l1}. Note that since we apply Scheffe's Lemma, we require the kernel $K_{\epsilon}$ to be non-negative. The bias analysis includes the following two cases. The proof of the proposition is in Appendix \ref{Appendix C}.

\begin{proposition}\label{bias L1 case}
\begin{enumerate}
\item
 Under Assumption \ref{assumptions on kernel}, suppose $P$ is continuous on $M$ with $\|P\|_\infty = P_{max}$  for some $P_{max}>0$,  then $\mathbb{E}K_{\epsilon}(x) \rightarrow P(x)$ for all $x \in M$ as $\epsilon \rightarrow 0$.
\item
 Under Assumption \ref{assumptions on kernel}, suppose $\alpha=d$ and $P$ is measurable on $M$ with $\|P\|_\infty = P_{max}$ for some $P_{max}>0$,  then $\mathbb{E}K_{\epsilon}(x) \rightarrow P(x)$ for  almost every $x \in M$ as $\epsilon \rightarrow 0$.
\end{enumerate}
\end{proposition}



Since we have
\begin{align}
&\int_{M}\left |\frac{1}{n}\sum_{i=1}^n K_{\epsilon}(\iota(x),\iota(x_i))-P(x)\right| dV(x)\nonumber \\
\leq &\,\int_{M} \left|\frac{1}{n}\sum_{i=1}^n K_{\epsilon}(\iota(x),\iota(x_i))-\mathbb{E}K_{\epsilon}(x)\right| dV(x) + \int_{M} |\mathbb{E}K_{\epsilon}(x)-P(x)| dV(x)\,, \nonumber
\end{align}
we can concatenate the previous two propositions and conclude the convergence of $K_{n,\epsilon}$ to $P$ in the $L^1$ sense.

\begin{theorem}\label{KDE L1 convergence}
Under Assumption \ref{assumptions on kernel}, suppose $\epsilon \rightarrow 0$ and $\frac{\log n}{n \epsilon^{2\alpha-d}} \rightarrow 0$ as $n \rightarrow \infty$. Then, a.s.
\begin{align}
\int_{M}|K_{n,\epsilon}(x)-P(x)|dV(x) \rightarrow 0
\end{align}
whenever one of the following conditions holds.
\begin{enumerate}
\item
$P$ is continuous on $M$ with $\|P\|_\infty = P_{max}$ for some $P_{max}>0$.
\item
$\alpha=d$ and $P$ is measurable on $M$ with $\|P\|_\infty = P_{max}$ for some $P_{max}>0$.
\end{enumerate}
\end{theorem}

In a specific case, by using the same argument as in Theorem \ref{KDE L1 convergence}, we can prove the convergence of the density estimation in $L^1$ sense for an isotropic nonnegative Lebesgue integrable kernel.

\begin{theorem}
Suppose $K(t):\mathbb{R}_{\geq 0} \rightarrow \mathbb{R}_{\geq 0}$ is a non-negative Lebesgue integrable function satisfying conditions (i), (iii) and (iv) in Assumption \ref{assumptions on kernel 1}, and $P$ is measurable on $M$ with $\|P\|_\infty = P_{max}$ for some $P_{max}>0$.  Suppose $\epsilon \rightarrow 0$ and $\frac{\log n}{n \epsilon^{d}} \rightarrow 0$ as $n \rightarrow \infty$. Then, for the kernel defined in \eqref{Definition Kn isotropic}, we have a.s.
\begin{align}
\int_{M}|K_{n}(x)-P(x)|dV(x) \rightarrow 0
\end{align}
\end{theorem}






\section*{Acknowledgements}
The authors would like to thank Professor Tyrus Berry and Professor Yen-Chi Chen about the helpful suggestions of the literatures. 
\



\bibliographystyle{plain}
