---
title: Fast Samplet Transform
url: https://www.emergentmind.com/topics/fast-samplet-transform
type: topic
---

# Fast Samplet Transform

Searching arXiv for the specified paper to ground the response and verify metadata.
arXiv search query: 2503.17487 Samplets Fast Samplet Transform Harbrecht Multerer
The Fast Samplet Transform (FST) is the change of basis associated with samplets, a wavelet-type system for scattered data introduced by Harbrecht and Multerer. Samplets are signed measures defined on sets of arbitrarily distributed data sites in possibly high dimension, and they transfer familiar wavelet concepts—fast basis transforms, data compression, operator compression, operator arithmetics, and sparse approximation—to non-grid data. In the same framework, samplet matrix compression facilitates rapid solution of scattered data interpolation problems, including kernels with nonlocal support, and sparsity constraints become meaningful in samplet coordinates [2503.17487].

## 1. Functional setting and definition

Let \(X=\{x_1,\dots,x_N\}\subset \Omega\subset \mathbb{R}^d\) be an arbitrary set of data sites. The construction is formulated in the space
$$
X=\operatorname{span}\{\delta_{x_i}: i=1,\dots,N\}\subset [C(\Omega)]',
$$
that is, the space of finitely supported signed measures. The inner product is prescribed by
$$
\langle \delta_{x_i},\delta_{x_j}\rangle_X=m_{ij},
$$
with the usual choice \(m_{ij}=\delta_{ij}\), so \((X,\langle\cdot,\cdot\rangle_X)\cong \mathbb{R}^N\) is a Hilbert space.

A multiresolution analysis in \(X\) is a nested sequence
$$
X_0\subset X_1\subset \cdots \subset X_J=X,
$$
with levelwise splittings
$$
X_{j+1}=X_j\oplus \mathcal{S}_j.
$$
Here \(X_j\) is the scaling space and \(\mathcal{S}_j\) is the detail space. A samplet basis consists of scaling distributions \(\{\phi_{0,k}\}\) spanning \(X_0\) together with samplets \(\sigma_{j,k}\) spanning each \(\mathcal{S}_j\), for \(j=0,\dots,J-1\).

The defining wavelet-type condition is polynomial cancellation. Each samplet satisfies \(q+1\) polynomial vanishing moments,
$$
\int_\Omega p(x)\,d\sigma_{j,k}(x)=0,\qquad \forall p\in \mathcal{P}_q,
$$
where \(\mathcal{P}_q\) denotes the polynomials of total degree at most \(q\), with
$$
\dim \mathcal{P}_q=m_q=\binom{q+d}{d}.
$$
This distinguishes samplets from arbitrary local signed measures: they are designed to annihilate low-order polynomial content while remaining localized on scattered point sets.

## 2. Multiresolution construction on cluster trees

The construction uses a balanced binary cluster tree \(\mathcal{T}\) of the data sites, for example by median cut. Each node \(\tau\subset X\) has two children, and refinement continues until the leaf size is \(\lesssim m_q\). If \(\operatorname{level}(\tau)=j\), then the node carries scaling distributions \(\Phi_j^\tau\) and samplets \(\Sigma_j^\tau\).

On each cluster \(\tau\) at level \(j\), the local two-scale relation is
$$
[\Phi_j^\tau,\Sigma_j^\tau]=\Phi_{j+1}^\tau\,Q_j^\tau,
$$
with
$$
Q_j^\tau=
\begin{bmatrix}
Q_{j,\Phi}^\tau\\
Q_{j,\Sigma}^\tau
\end{bmatrix}
\in \mathbb{R}^{n_{j+1}^\tau\times n_{j+1}^\tau},
\qquad n_{j+1}^\tau=|\Phi_{j+1}^\tau|.
$$
The role of \(Q_j^\tau\) is to separate scaling and detail components while enforcing orthonormality and vanishing moments.

The local moment matrix is
$$
M_{j+1}^\tau=
\big[\langle x^\alpha,\phi_{j+1,\ell}^\tau\rangle_\Omega\big]_{|\alpha|\le q,\ \ell=1,\dots,n_{j+1}^\tau}.
$$
In practice one computes
$$
(M_{j+1}^\tau)^T=Q_j^\tau R_j
$$
by Gram–Schmidt or Householder QR. The first \(m_q\) columns of \(Q_{j,\Phi}^\tau\) define the scaling part, and the remaining \(n_{j+1}^\tau-m_q\) columns of \(Q_{j,\Sigma}^\tau\) define the samplets.

Several structural properties follow directly from this construction. The family \(\{\phi_{0,k}\}\cup\{\sigma_{j,k}\}\) is an orthonormal basis of \(X\). The number of samplets at level \(j\) is \(\simeq 2^j\), with total cardinality \(N\). Each \(\sigma_{j,k}\) is supported in one cluster \(\tau\) and has \(q+1\) vanishing moments. A common misconception is to treat samplets as ordinary wavelet functions on a regular mesh; in this framework they are signed measures supported on clustered scattered data sites.

## 3. Forward and inverse Fast Samplet Transform

The FST is the matrix \(T\in \mathbb{R}^{N\times N}\) defined by
$$
f^\Sigma=T\,f^\Delta,
$$
where
$$
f_i^\Delta=f(x_i),\qquad f_{j,k}^\Sigma=\int_\Omega f\,d\sigma_{j,k}.
$$
The transform matrix is block-sparse, with a “fish-bone” structure, and is applied through local two-scale relations rather than dense matrix multiplication.

At cluster \(\tau\) and level \(j\), the decomposition identity is
$$
[(f,\phi_j^\tau)_\Omega,\ (f,\sigma_j^\tau)_\Omega]
=
(f,\Phi_{j+1}^\tau)_\Omega\,[Q_{j,\Phi}^\tau,Q_{j,\Sigma}^\tau].
$$
The forward transform starts from the finest-scale values
$$
f_J^\Phi=f^\Delta=[f(x_1),\dots,f(x_N)]^T,
$$
and processes levels \(j=J-1,J-2,\dots,0\). On each cluster \(\tau\), the restricted vector \(f_{j+1}^\Phi|_\tau\) is multiplied by
$$
[Q_{j,\Phi}^\tau,Q_{j,\Sigma}^\tau]^T
$$
to produce a new scaling part \(f_j^\Phi\) and a detail part \(f_j^\Sigma\). The detail coefficients are collected into \(f^\Sigma\), and the remaining \(f_0^\Phi\) is stored as the coarsest-scale component.

Reconstruction reverses the process. At each level, the coarse scaling data and detail data on a cluster are combined using
$$
[Q_{j,\Phi}^\tau,Q_{j,\Sigma}^\tau]
$$
to recover the next finer-scale vector, culminating in \(f_J^\Phi=f^\Delta\).

In block form, the levelwise analysis and synthesis filters are
$$
A_j=\bigoplus_{\tau\text{ at level }j}[Q_{j,\Phi}^\tau,Q_{j,\Sigma}^\tau]^T,
\qquad
S_j=\bigoplus_{\tau\text{ at level }j}[Q_{j,\Phi}^\tau,Q_{j,\Sigma}^\tau].
$$
Hence
$$
f^\Sigma=A_0A_1\cdots A_{J-1}f^\Delta,
\qquad
f^\Delta=S_{J-1}\cdots S_1S_0\,f^\Sigma.
$$
This formulation makes explicit that the FST is a hierarchical exact basis change before any compression is applied.

## 4. Computational complexity and coefficient compression

The cluster tree can be built in \(O(N\log N)\) time by median-split on bounding boxes [2503.17487]. For the transform itself, the depth satisfies \(J=O(\log N)\). At each cluster \(\tau\), the local work involves a small block of size governed by \(m_q\), and summing over all levels yields
$$
O(m_q^2N).
$$
Since
$$
m_q=\binom{q+d}{d}=O(d^q),
$$
both analysis and synthesis are \(O(N)\) for fixed \(q\) and \(d\). The same linear bound applies to coefficient thresholding, and if only nonzero coefficients are processed after thresholding, the cost reduces proportionally to the number of nonzeros.

Coefficient compression is implemented by hard thresholding. For a threshold \(w>0\),
$$
HT_w(v)_i=
\begin{cases}
v_i,& |v_i|\ge w,\\
0,& \text{otherwise}.
\end{cases}
$$
Compression is therefore obtained by replacing \(f^\Sigma\) with \(HT_w(f^\Sigma)\) before synthesis. This separation between exact multiresolution analysis and optional nonlinear truncation parallels standard wavelet workflows, but here it is realized on arbitrarily distributed sites rather than on a structured grid.

## 5. Tuning parameters and operational tradeoffs

The principal tuning parameter is \(q+1\), the number of vanishing moments. Larger \(q\) yields faster decay of \(f^\Sigma\) for smooth \(f\), but also enlarges the local block size \(m_q\). Typical values are \(q+1=3\ldots 5\). This is the main approximation-versus-local-cost tradeoff built into the construction.

Leaf size should satisfy \(\gtrsim m_q\); examples given are \(\mathrm{leaf\_size}=2m_q\) or a constant in the range \(20\ldots 50\). Thresholding uses the parameter \(w\), and a common practical choice is
$$
w=\epsilon\cdot \|f^\Sigma\|_2,
$$
with \(\epsilon=10^{-3}\ldots 10^{-5}\). For kernel-matrix compression, an admissibility parameter \(\eta>0\) is introduced. Tree balance need not be exact: the cost bounds continue to hold for approximately balanced trees.

These parameters control distinct aspects of the method. The moment order controls cancellation and coefficient decay, leaf size controls the local algebra, the threshold controls nonlinear compression in samplet coordinates, and \(\eta\) controls admissibility in kernel-matrix compression. Within the stated ranges, the chapter reports effective behavior on data sets from \(N\approx 10^6\) up to \(N=1.2\) million.

## 6. Reported numerical behavior

The numerical experiments span coefficient compression, kernel-matrix compression, scattered-data interpolation, and sparse basis pursuit [2503.17487].

| Problem | Setup | Reported outcome |
|---|---|---|
| Hard-thresholding of global temperature (ERA5) | Monthly data, \(N\approx 10^6\) per month, \(q+1=4\), \(w=10^{-k}\|f\|_2\), \(k=2,3,4,5\) | Nonzeros \(\approx 4,87,1\,454,15\,688\); space saving \(99.9996\%,99.9916\%,99.86\%,98.49\%\); relative \(\ell^2\)-errors \(2.4\%,1.2\%,0.39\%,0.11\%\) |
| Kernel-matrix compression | Exponential kernel \(k(r)=e^{-r/\ell}\) on unit square, \(\ell=0.1\), \(q+1=4\), \(\eta=1\) | \(\|K^\Sigma-K_\epsilon^\Sigma\|_F/\|K^\Sigma\|_F=O((\eta\rho)^{-2(q+1)})\) with \(O(N\log N)\) nonzeros; after nested-dissection reordering, an extremely sparse Cholesky factor with near-linear fill |
| Surface reconstruction (“Laokoon” signed-distance) | \(N=750\,000\) points, \(\ell=0.01\), \(q+1=4\), \(\mu=10^{-8}\) Tikhonov | Build and compress \(K\) in \(O(N\log N)\); solve via CG in \(\approx 100\,\mathrm{s}\); evaluate on \(125\) million grid points via FMM in \(\approx 50\,\mathrm{s}\) |
| Sparse basis pursuit on space–time temperature | Two kernels \(k_{3/2}\otimes k_{\mathrm{per}}\) and \(k_{1/2}\), \(N=1.2\) million adaptively subsampled points, \(q+1=4\), \(w_i=10^{-6}\) uniform penalty | Semi-smooth Newton in compressed format; \(\|\beta_1\|_0=671\), \(\|\beta_2\|_0=5\,883\), in-sample \(\|h-K\alpha\|_2/\|h\|_2\approx 8\times 10^{-3}\), max pointwise relative error \(<9.8\%\) |

For the ERA5 thresholding experiment, the reported pointwise maximum error is up to \(21\%\) for the coarsest threshold and \(\lesssim 2\%\) for the finest. For the surface-reconstruction experiment, the reported outcome is visually high-quality, smooth level-sets together with interactive run-times on large \(N\). Across the examples, the recurring pattern is that sparsity in samplet coordinates is operationally useful for both vectors and kernel matrices.

## 7. Conceptual position within scattered-data analysis

The chapter places samplets within wavelet analysis for scattered data rather than within classical grid-based wavelet theory. The essential shift is from basis functions tied to regular meshes to signed measures supported on local clusters of arbitrary data sites. This permits multiresolution structure, vanishing moments, fast transforms, and threshold-based compression without requiring geometric regularity of the sampling set.

The reported uses are correspondingly broad: fast basis transform, data compression, operator compression, operator arithmetics, rapid interpolation for kernels with nonlocal support, and sparse approximation under explicit sparsity constraints in samplet coordinates. The further-reading list associated with the chapter includes work by Beylkin, Coifman, and Rokhlin; Tausch and White; and later work by Harbrecht, Multerer, Schenk, and Schwab. This suggests a lineage connecting classical wavelet compression, multiresolution methods on unstructured data, and kernel-algebra techniques.

Two points delimit the scope of the method clearly. First, samplets are not merely compressed point values; they are basis elements with prescribed \(q+1\) polynomial vanishing moments and cluster support. Second, the FST is not only a compression heuristic; it is an \(O(N)\) analysis/synthesis mechanism for an orthonormal basis, with hard thresholding added as a separate nonlinear step. In that sense, the Fast Samplet Transform is best understood as the multiresolution computational core of the samplet framework for scattered data.

Source: https://www.emergentmind.com/topics/fast-samplet-transform