Kernel Interpolation in TP-RKHS

• Kernel interpolation in TP-RKHSs is a framework for high-dimensional function approximation using tensor-product kernels that ensure strict positive definiteness.
• The approach leverages structured grids and Kronecker product factorization to reduce computational complexity in interpolating multivariate functions.
• Optimized sparse grid techniques coupled with hybrid Sobolev norms eliminate logarithmic penalties, recovering univariate convergence rates in high dimensions.

Kernel interpolation in tensor product reproducing kernel Hilbert spaces (TP-RKHSs) addresses the problem of high-dimensional function approximation and interpolation by leveraging kernel-based methods structured over tensor-product domains. This approach exploits fundamental properties of RKHSs, especially when constructed as tensor products, to achieve scalable interpolation algorithms, error bounds reflecting univariate convergence rates, and, under refined sparse grid constructions, the possibility of entirely avoiding logarithmic penalty factors classically associated with high dimension.

1. Foundational Structures: Tensor Product RKHSs and Product Kernels

Let $H^p(\mathbb T)$ denote the $p$-smooth periodic Sobolev space with norm induced by an inner product encoding Fourier coefficients, and $k_p(x,y)$ the associated reproducing kernel. For domains $\Omega_j \subset \mathbb R$ and univariate kernels $k_j$, the construction of a tensor product kernel proceeds via

$k(x,y) = \prod_{j=1}^d k_j(x_j, y_j)$

on $\Omega = \Omega_1 \times \cdots \times \Omega_d$.

The native space of this product kernel is isometrically isomorphic to the Hilbert tensor product $\mathcal H_1 \otimes \cdots \otimes \mathcal H_d$, with evaluation functional $$f(x) = \langle f, k(x,\,\cdot\,) \rangle_{\mathcal H_{k,\Omega}}$$ and inner product

$$\langle \varphi(f_1, \dots, f_d), \varphi(g_1, \dots, g_d) \rangle = \prod_{j=1}^d \langle f_j, g_j \rangle_{\mathcal H_j}$$

(Albrecht et al., 2023). This structure guarantees strict positive definiteness of $k$ provided each $k_j$ is strictly positive definite, which is essential for uniqueness in interpolation problems.

For $d$-fold periodic domains, the mixed Sobolev space $H^p_{mix}(\mathbb T^d)$ coincides with the $d$-fold tensor product $\bigotimes_{i=1}^d H^p(\mathbb T)$, with the kernel $$K(\mathbf x, \mathbf y) = \prod_{i=1}^d k_p(x_i, y_i)$$ (Griebel et al., 14 Dec 2025).

2. Kernel Interpolation Algorithms in Tensor Product Spaces

Given a set of nodes $X = \{x_i\}_{i=1}^n \subset \Omega$ and corresponding data values $f_i$, the interpolant $s$ in the span of $\{k(\cdot, x_i)\}$ is determined by the system

$K\alpha = f$

with $K_{ij} = k(x_i, x_j)$. For tensor-product (grid) nodes $X_1 \times \cdots \times X_d$, the Gram matrix factorizes as a Kronecker product $K = K^{(1)} \otimes \cdots \otimes K^{(d)}$, significantly reducing computational complexity:

• 1D Cholesky: $O(n_j^3)$,
• $d$-dimensional grid: $\sum_j O(n_j^3)$ vs $O(N^3)$ for naive matrices, where $N = \prod_j n_j$, enabling large-scale interpolation schemes via Kronecker and FFT/block-Krylov structures (Albrecht et al., 2023, Griebel et al., 14 Dec 2025).

For scattered data, the Newton basis can be constructed univariately using Cholesky factorizations, and then multivariately as tensor products, directly reflecting the basis structure of the tensor-product RKHS (Albrecht et al., 2023).

3. Sparse Grids and Hybrid Regularity: Dimensionality Reduction

Classical tensor grid methods for interpolation suffer from the curse of dimensionality: full grids have exponentially many nodes, and classical sparse grids (Smolyak-type, $\|\mathbf j\|_1 \le J$) yield a log-factor penalty in convergence rates:

$\|u - u_J\| \lesssim N^{-r} (\log N)^{(d-1)(r+1)}$

in $H^t_{mix}$ (Griebel et al., 14 Dec 2025).

By contrast, choosing an optimized sparse grid index set

$$\mathcal I_J^{\lambda} = \{\mathbf j \in \mathbb N_0^d: \|\mathbf j\|_1 - \lambda \|\mathbf j\|_\infty \le J(1-\lambda)\}$$

with any $\lambda \in (0, (s_1-s_2)/(t_2-t_1))$, and measuring error in Sobolev spaces of hybrid regularity $H_{iso\text{-}mix}^{s,t}$, the logarithmic penalty disappears:

$$\|u - \widehat Q_J^\lambda u\|_{H_{iso\text{-}mix}^{s_1,t_1}} \le C N^{-r} \|u\|_{H_{iso\text{-}mix}^{s_2,t_2}}$$

for $N \sim 2^J$, with $r = (t_2-t_1) - (s_1-s_2) > 0$ and no dependence on $\log N$ or $d$ beyond the constant (Griebel et al., 14 Dec 2025). This recovers univariate convergence rates in high-dimensional settings.

The hybrid regularity space $H_{iso\text{-}mix}^{s,t}(\mathbb T^d)$ consists of functions whose derivatives in each coordinate direction enjoy extra $s$-regularity beyond their mixed-regularity $t$, yielding embeddings that are tight enough for the optimized sparse grid analysis.

4. Theoretical Guarantees and Error Estimates

In the optimized setting,

• Recovery rate: Univariate-style $N^{-r}$ convergence for functions in a Sobolev space of hybrid regularity.
• Complexity: $O(N \log N)$ total work for fixed $d$, when using Kronecker/FFT methods for each grid component.
• Best-possible grids: For $\lambda > 0$ (i.e., non-classical sparse grids), the number of degrees of freedom remains $O(2^J)$, not $O(2^{Jd})$.
• No logarithmic loss: Provided interpolation is measured in $H_{iso\text{-}mix}^{s_1,t_1}$ with extra smoothness as in the analysis.

Comparative results show classical approaches (Smolyak, Korobov/isotropic Sobolev spaces) are subject to unavoidable $(\log N)^{d-1}$ factors unless highly specialized periodic kernels and custom combination weights are used (Griebel et al., 14 Dec 2025).

5. Multilinear Spectral Penalization in TP-RKHS Interpolation

For general M-way TP-RKHSs $\mathcal H_1 \otimes \cdots \otimes \mathcal H_M$, interpolation problems can be formulated with additional structure-promoting regularization:

• Given training points $\{x^{(i)}\}_{i=1}^n$ and labels $y^{(i)}$, functions $f$ in TP-RKHS admit representations via coefficient tensors $A$:

$$f(x) = \sum_{i_1, \ldots, i_M} a_{i_1 \ldots i_M} \prod_{m=1}^M k_m(x^m, x^{(i_m),m})$$

Constrained exact interpolation is encoded as

$$A \times_1 K^{(1)} \times_2 \cdots \times_M K^{(M)} = Y$$

for Gram matrices $K^{(m)}$, vector $Y$ of targets (Signoretto et al., 2013).

Spectral penalties, particularly the sum of nuclear norms of $A$'s mode-wise unfoldings,

$$\Omega_{MSP}(A) = \sum_{m=1}^M \gamma_m \|A_{(m)}\|_*,$$

encourage low-multilinear-rank solutions, implicitly regularizing the complexity of the interpolant in each direction. The corresponding optimization is a convex program admitting a block-coordinate SVD-proximal solution or Tucker-decomposition-based alternating minimization (Signoretto et al., 2013).

A plausible implication is that, for datasets with underlying low-rank tensor structure, such penalized TP-RKHS interpolation schemes deliver both statistical and algorithmic advantages.

6. Practical Implementation and Algorithmic Considerations

Implementation for kernel interpolation in TP-RKHS with optimized sparse grids entails:

• Selection of a univariate kernel $k_p$ and determination of suitable regularity indices $(s_1, t_1)$ and $(s_2, t_2)$;
• Choosing $\lambda$ and $J$ so $N \sim 2^J$ matches computational or accuracy budget;
• Constructing $\mathcal I_J^\lambda$ and, for each $\mathbf j$, computing $P_{\mathbf j} u$ via Kronecker/FFT/block-Krylov methods;
• Aggregating results with Boolean weights $c_{\mathbf j}$:

$$\widehat u_J = \sum_{\mathbf j \in \mathcal I_J^\lambda} c_{\mathbf j} P_{\mathbf j} u$$

achieving error bounds without logarithmic factors (Griebel et al., 14 Dec 2025).

For scattered data or non-grid settings, explicit use of the Newton basis and the Kronecker structure can yield order-of-magnitude computational benefits. Product kernels enable blending properties (stability, sharpness) across dimensional axes, as supported by empirical timing and error studies on two-dimensional test functions (Albrecht et al., 2023).

In cases where multilinear spectral penalties are employed, the core computation involves iterated SVD and dual ascent steps, taking particular advantage of tensor factorizations and orthogonal projections (Signoretto et al., 2013).

7. High-Dimensional Interpolation: Impact and Limitations

The developments in kernel interpolation over TP-RKHSs, especially when combined with optimized sparse grid constructions and hybrid Sobolev regularity, permit univariate-type convergence rates in genuinely high-dimensional approximation tasks. This mitigates the curse of dimensionality without recourse to exotic kernel choices or problem-dependent weights for the combination formula.

The absence of log-factors in the overall error-to-complexity relationship is strictly contingent on measuring error in hybrid regularity norms and employing non-classical (optimized) sparse index sets. In strictly isotropic or pure mixed regularity settings, or with classical sparse grids, logarithmic losses remain.

A plausible implication is that for high-dimensional regression, scattered data approximation, and learning problems with underlying hybrid regularity, the combination of TP-RKHS methodology, optimized sparse grids, and possibly multilinear spectral regularization forms a theoretically optimal and computationally tractable paradigm (Griebel et al., 14 Dec 2025, Albrecht et al., 2023, Signoretto et al., 2013).