1 Reproducing kernel Hilbert space

Let $X\subset {ℝ}^d$ be a non-empty set. A function $K:X\times X\to ℝ$ is called a symmetric positive definite (SPD) kernel if

$K(x,y)=K(y,x)\mathrm{and}{\displaystyle \sum _{i,j=1}^n}{c}_i{c}_{j}K(x_i,x_j)\ge 0$ (1)

for every finite set of points $\{x_1,\mathrm{\dots },x_n\}\subset X$ and every choice of real coefficients $\{c_1,\mathrm{\dots },c_n\}\subset ℝ$.

The classical Moore-Aronszajn theorem asserts that for every such kernel $K$ there exists a unique Hilbert space $(\mathscr{H},{⟨\cdot ,\cdot ⟩}_{\mathscr{H}})$ of real-valued functions on $X$ with the following properties:

1. For each $x\in X$ the function $K(\cdot ,x)$ belongs to $\mathscr{H}$;

2. For every $f\in \mathscr{H}$ and every $x\in X$,

   $f(x)={⟨f,K(\cdot ,x)⟩}_{\mathscr{H}}.$ (2)

The Hilbert space $\mathscr{H}$ is called the reproducing kernel Hilbert space (RKHS) associated with the kernel $K$. The function $K$ itself is the reproducing kernel of $\mathscr{H}$. Property (2) characterizes the evaluation functional as a continuous linear functional on $\mathscr{H}$, which is usually called the reproducing property.

2 RKHS with Matérn kernel and Sobolev space

A particularly important class of RKHSs arises from shift-invariant kernels of the form

$$K(x,y)=\mathrm{\Phi }(x-y),\mathrm{\Phi }:{ℝ}^d\to ℝ,$$

where $\mathrm{\Phi }$ is an SPD kernel function. The next result describes the structure of the RKHS associated with such a kernel; see [4, Theorem 10.12].

Theorem 1 (RKHS of a shift-invariant kernel)

Let $K(x,y)=\mathrm{\Phi }(x-y)$ be a shift-invariant kernel on ${ℝ}^d$ with $\mathrm{\Phi }\in C({ℝ}^d)\cap L^1({ℝ}^d)$. Then the associated RKHS is

$\mathscr{H}=\{f\in L^2({ℝ}^d)\cap C({ℝ}^d):ℱ{[\mathrm{\Phi }]}^{-\frac{1}{2}}ℱ[f]\in L^2({ℝ}^d)\},$ (3)

where $ℱ[\mathrm{\Phi }]$ is the Fourier transform of $\mathrm{\Phi }$. The inner product on $\mathscr{H}$ is

$${⟨f,g⟩}_{\mathscr{H}}=\frac{1}{{(2\pi )}^{d/2}}{\int }_{{ℝ}^d}\frac{ℱ[f](\omega )\overline{ℱ[g](\omega )}}{ℱ[\mathrm{\Phi }](\omega )}d\omega .$$

Consider the Matérn kernel

$K_{\nu ,l}(x,y)={\mathrm{\Phi }}_{\nu ,l}(x-y),{\mathrm{\Phi }}_{\nu ,l}(x)={\displaystyle \frac{2^{1-\nu }}{\mathrm{\Gamma }(\nu )}}{\left(\sqrt{2\nu }{\displaystyle \frac{\Vert x\Vert }{l}}\right)}^{\nu }{B}_{\nu }\left(\sqrt{2\nu }{\displaystyle \frac{\Vert x\Vert }{l}}\right),$ (4)

where $\nu ,l>0$, and $B_{\nu }$ is the modified Bessel function of the second kind. The Fourier transform of ${\mathrm{\Phi }}_{\nu ,l}$ is

$ℱ[{\mathrm{\Phi }}_{\nu ,l}](\omega )=C_{\nu ,l,d}{\left({\displaystyle \frac{2\nu }{l^2}}+4{\pi }^2{\Vert \omega \Vert }^2\right)}^{-(\nu +\frac{d}{2})},$ (5)

where $C_{\nu ,l,d}$ is a constant only depending on $(\nu ,l,d)$; see [5, §4.2]. Thus, the condition in (3) becomes ${\left(\frac{2\nu }{l^2}+4{\pi }^2{\Vert \omega \Vert }^2\right)}^{(\nu +\frac{d}{2})/2}\widehat{f}(\omega )\in L^2({ℝ}^d)$. Recall that for any $s\in ℝ$, the fractional Sobolev space on ${ℝ}^d$ is

$H^s({ℝ}^d)=\{f\in L^2({ℝ}^d):{(1+{\Vert \omega \Vert }^2)}^{\frac{s}{2}}\widehat{f}(\omega )\in L^2({ℝ}^d)\},$ (6)

where $\widehat{f}(\omega )={(2\pi )}^{-d/2}{\int }_{{ℝ}^d}f(x)e^{-i⟨\omega ,x⟩}dx$ is the Fourier transform of $f$. The Sobolev embedding theorem states that $H^s({ℝ}^d)\subset C({ℝ}^d)$ if $s>\frac{d}{2}$. Note that $c_1(1+{\Vert \omega \Vert }^2)\le \left(\frac{2\nu }{l^2}+4{\pi }^2{\Vert \omega \Vert }^2\right)\le c_2(1+{\Vert \omega \Vert }^2)$ with some $c_1,c_2>0$. Therefore, combining (3), (5) and (6), we conclude that: the RKHS with Matérn kernel $K(x,y)={\mathrm{\Phi }}_{\nu ,l}(x-y)$ is norm-equivalent to the Sobolev space $H^{\nu +\frac{d}{2}}({ℝ}^d)$.

For a domain $\mathrm{\Omega }\in {ℝ}^d$ and $s\ge 0$, the Sobolev space $H^s(\mathrm{\Omega })$ is the set of restrictions of functions from $H^s({ℝ}^d)$ to $\mathrm{\Omega }$ equipped with the norm

$\parallel f{\parallel }_{H^s(\mathrm{\Omega })}:=inf\{\parallel g{\parallel }_{H^s({ℝ}^d)}:g\in H^s({ℝ}^d),g{|}_{\mathrm{\Omega }}=f\}.$

The RKHS with kernel $K_{\nu ,l}$ can also be defined on $\mathrm{\Omega }$. If $\mathrm{\Omega }$ is a Lipschitz domain, then $\mathscr{H}$ with Matérn kernel is norm-equivalent to the Sobolev space $H^{\nu +\frac{d}{2}}(\mathrm{\Omega })$; see [4, §10.7].

3 Gaussian random field and covariance operator

The following result states that a Gaussian random field (GRF) with Matérn kernel can be given as the solution of a stochastic partial differential equation (SPDE); see [2].

Theorem 2

Let $GP(0,K_{\nu ,l})$ be a centered GRF on a Lipschitz domain $\mathrm{\Omega }\subset {ℝ}^d$ such that the covariance function is the Matérn kernel $K_{\nu ,l}$, and $\mathrm{\Delta }$ be the Laplacian with Dirichlet boundary condition. Then $u\sim GP(0,K_{\nu ,l})$ is the unique solution of the fractional elliptic SPDE

${\left({\displaystyle \frac{2\nu }{l^2}}-\mathrm{\Delta }\right)}^{\frac{\gamma }{2}}u=W,\gamma =\nu +{\displaystyle \frac{d}{2}},$ (7)

where $W$ is a spatial Gaussian white noise on $L^2(\mathrm{\Omega })$.

The above SPDE is defined in the distributional sense. Note that $W$ is an isometry from $L^2(\mathrm{\Omega })$ to a centered Gaussian space and write $⟨W,\varphi ⟩:=W(\varphi )$. Since ${(2\nu /l^2-\mathrm{\Delta })}^{-\frac{\gamma }{2}}$ is a Hilbert-Schmidt operator with $\gamma >\frac{d}{2}$, this implies that the random element $u$ lies in $L^2(\mathrm{\Omega })$ almost surely, and the GRF $GP(0,K_{\nu ,l})$ induces a Gaussian measure on $L^2(\mathrm{\Omega })$. The Cameron-Martin space of the induced Gaussian measure is the RKHS ${\mathscr{H}}_{K_{\nu ,l}}$, which is norm-equivalent to $H^{\nu +\frac{d}{2}}(\mathrm{\Omega })$ and continuously embedded into $L^2(\mathrm{\Omega })$. For any $\varphi \in L^2(\mathrm{\Omega })$, the SPDE implies that

$⟨W,\varphi ⟩={⟨{(2\nu /l^2-\mathrm{\Delta })}^{\frac{\gamma }{2}}u,\varphi ⟩}_{L^2(\mathrm{\Omega })}={⟨u,{(2\nu /l^2-\mathrm{\Delta })}^{\frac{\gamma }{2}}\varphi ⟩}_{L^2(\mathrm{\Omega })},$

hence $⟨W,{(2\nu /l^2-\mathrm{\Delta })}^{-\frac{\gamma }{2}}\psi ⟩={⟨u,\psi ⟩}_{L^2(\mathrm{\Omega })}$ for any $\psi \in L^2(\mathrm{\Omega })$ since $(2\nu /l^2-\mathrm{\Delta })$ is invertible on $L^2(\mathrm{\Omega })$. For any $\varphi ,\psi \in L^2(\mathrm{\Omega })$, it follows that

$𝔼[{⟨u,\varphi ⟩}_{L^2(\mathrm{\Omega })}{⟨u,\psi ⟩}_{L^2(\mathrm{\Omega })}]$ $=𝔼[⟨W,{(2\nu /l^2-\mathrm{\Delta })}^{-\frac{\gamma }{2}}\varphi ⟩⟨W,{(2\nu /l^2-\mathrm{\Delta })}^{-\frac{\gamma }{2}}\psi ⟩]$

$={⟨{(2\nu /l^2-\mathrm{\Delta })}^{-\frac{\gamma }{2}}\varphi ,{(2\nu /l^2-\mathrm{\Delta })}^{-\frac{\gamma }{2}}\psi ⟩}_{L^2(\mathrm{\Omega })}$

$={⟨{(2\nu /l^2-\mathrm{\Delta })}^{-\gamma }\varphi ,\psi ⟩}_{L^2(\mathrm{\Omega })}.$

Note that by the Fernique’s theorem. This implies that $B_u$ define by $B_u(\varphi )={⟨u,\varphi ⟩}_{L^2(\mathrm{\Omega })}$ is a regular zero-mean generalized Gaussian field on $L^2(\mathrm{\Omega })$, and its covariance operator is ${(2\nu /l^2-\mathrm{\Delta })}^{-\gamma }$. Therefore, the covariance operator can be written as the kernel integral operator

$(T_{K}f)(x)={\displaystyle {\int }_{\mathrm{\Omega }}}{K}_{\nu ,l}(x,y)f(y)dy$ (8)

defined on $L^2(\mathrm{\Omega })$; see [3, Exercise 3.2.14]. That is, it holds $T_K={(2\nu /l^2-\mathrm{\Delta })}^{-\gamma }$.

A byproduct of the above result is that we can get the decay rate of the eigenvalues of the trace-class operator $T_K$ on $L^2(\mathrm{\Omega })$ or $H^{\nu +\frac{d}{2}}(\mathrm{\Omega })$ (the eigenvalues are the same on these two spaces). For a bounded domain $\mathrm{\Omega }\subset {ℝ}^d$, recall that the Laplacian $-\mathrm{\Delta }$ with Dirichlet boundary condition has eigenvalues ${\{{\lambda }_j\}}_j$ (arranged in ascending order) that increase obeying the Weyl’s law ${\lambda }_j\asymp j^{2/d}$; see [1, §6.4]. Therefore, the eigenvalues of $T_K$ are the eigenvalues of ${(2\nu /l^2-\mathrm{\Delta })}^{-\gamma }$, which decays with the rate

${\mu }_j={(2\nu /l+{\lambda }_j)}^{-\gamma }\asymp j^{-2\gamma /d}=j^{-\left(1+\frac{2\nu }{d}\right)},j\to \mathrm{\infty }.$ (9)

Therefore, the eigenvalues of $T_K$ with Matérn kernel has a polynomial decay rate.