Source-linked AI summary
Stochastic trace estimation for positive trace-class operators
Zvonimir Bujanović, Luka Grubišić, Daniel Kressner, Hrvoje Olić
TL;DR
Implicit trace estimation needs efficient methods for operators accessible only through products, but standard and variance-reduced matrix estimators do not directly address the infinite-dimensional setting. The paper develops operator analogues with target-operator-dependent Gaussian covariances and finite-dimensional truncations, proving unbiasedness, high-probability guarantees, and convergence results. The truncated methods achieve accuracy comparable to ContHutch++ in numerical experiments while using lower-degree representations and smaller chebfun discretizations.
Problem
Implicit trace estimation for operators accessible only through operator-vector products lacks corresponding infinite-dimensional Girard–Hutchinson and Hutch++ methods with the matrix estimators’ guarantees.
Method
The paper constructs idealized estimators for positive trace-class operators using Gaussian random elements with covariance determined by the target operator, then truncates sampling to finite-dimensional subspaces.
Results
Truncated infHutch++ achieves similar or better stochastic accuracy than ContHutch++ with substantially smaller internal chebfun discretizations, while idealized infHutch++ has O(ε^-1) sample complexity.
Takeaways & Limitations
Operator-dependent covariance removes the smoothing bias associated with ContHutch++’s externally prescribed covariance, and truncated estimators provide a practical route toward the idealized methods.
Takeaways & Limitations
The adaptive truncation stopping rule is an empirical heuristic and does not provide a high-probability bound for the full trace tail without additional decay assumptions.
Abstract
from arXiv · showhide
Implicit trace estimation aims to approximate the trace of a matrix or linear operator accessible only through matrix-vector or operator-vector products. In the matrix setting, the Girard-Hutchinson estimator typically requires $\mathcal{O}(\varepsilon^{-2})$ products to achieve accuracy $\varepsilon$, while the variance-reduced Hutch++ algorithm reduces this sample complexity to $\mathcal{O}(\varepsilon^{-1})$ for positive semidefinite matrices. We develop infinite-dimensional analogues of these estimators for positive trace-class operators on separable Hilbert spaces. The idealized estimators use Gaussian random elements whose covariance is determined by the target operator, leading to unbiased operator versions of Girard-Hutchinson and Hutch++. We prove high-probability error bounds analogous to the finite-dimensional matrix results; in particular, idealized infHutch++ achieves $\mathcal{O}(\varepsilon^{-1})$ sample complexity. For practical computation, we introduce truncated implementations that restrict the random samples to finite-dimensional subspaces; for fixed sample budget, we show that truncated infHutch++ converges in distribution to its idealized counterpart as the truncation dimension tends to infinity. Numerical experiments with integral operators, density-of-states approximations, and spectral filtering for a radial Dirac operator show that these truncated estimators can achieve accuracy comparable to the ContHutch++ algorithm by Zvonek, Horning & Townsend while using lower-degree function representations and smaller internal discretizations in chebfun.
1 Introduction
The paper extends implicit stochastic trace estimation from matrices to positive trace-class operators on separable Hilbert spaces. It develops unbiased operator estimators, practical truncations, and analyses their theoretical and computational trade-offs.
- Implicit trace estimation approximates tr(A) using matrix-vector products and supports applications including spectral densities, log-determinants, and network measures.
- The Girard–Hutchinson estimator has Monte Carlo sample complexity O(ε^-2), whereas Hutch++ reduces this to O(ε^-1) for positive semidefinite matrices.
- The paper develops idealized Girard–Hutchinson and Hutch++ analogues for positive trace-class operators on separable Hilbert spaces.
- Idealized operator estimators use operator-dependent Gaussian distributions, while basis-truncated implementations address their lack of direct practical availability.
- Unlike ContHutch++, the proposed covariance operators are derived from A, avoiding externally tuned smoothing bias and length-scale regularity choices.
- The paper establishes infinite-dimensional theory, analyzes truncation consistency, and presents numerical experiments comparing the approach with ContHutch++.
2 Preliminaries
The preliminaries define the operator setting and review matrix trace estimation. They introduce trace-class operators, quasimatrices, Girard–Hutchinson, and Hutch++ as foundations for the infinite-dimensional constructions.
- Trace-class operators on a Hilbert space: The setting is a real, separable, infinite-dimensional Hilbert space containing positive trace-class operators.
- Trace-class operators on a Hilbert space: A positive trace-class operator is bounded, self-adjoint, positive, and has finite trace independent of the chosen orthonormal basis.
- Trace-class operators on a Hilbert space: Integral operators with continuous, symmetric, positive semidefinite kernels provide a typical example of the operator class.
- Quasimatrices: Quasimatrices represent finitely many Hilbert-space vectors as columns and provide an infinite-dimensional analogue of tall-and-skinny matrices.
- Matrix trace estimators: For isotropic random vectors, the Girard–Hutchinson estimator is unbiased, with nonasymptotic Gaussian tail bounds controlling its error.
- Matrix trace estimators: Hutch++ combines randomized low-rank approximation with stochastic trace estimation, reducing variance and achieving O(ε^-1) sample complexity for positive semidefinite matrices.
3 Operator trace estimation
The paper constructs unbiased Gaussian trace estimators for positive trace-class operators and extends Hutch++ through randomized low-rank approximation. Truncated implementations converge to the idealized infHutch++ estimator as the truncation dimension grows, while adaptive truncation remains heuristic.
- Idealized infHutch: Gaussian random elements with covariance A define an operator analogue of Girard–Hutchinson, called infHutch, whose expectation equals tr(A).The construction uses the spectral decomposition of A and establishes the Gaussian distribution N(0,A) in the Hilbert space.
- Idealized infHutch: For fixed failure probability, infHutch requires O(ε^-2) samples to achieve relative accuracy ε.The bound follows from concentration of ||x||^2 − tr(A) and trace-class norm inequalities.
- The infHutch++ algorithm: infHutch++ combines a randomized low-rank approximation with unbiased residual trace estimation using Gaussian samples from N(0,A).The range-finding samples come from N(0,A^2), while the residual is estimated by projecting samples from N(0,A) outside the computed subspace.
- The infHutch++ algorithm: The idealized infHutch++ estimator reduces the sample complexity from O(ε^-2) to O(ε^-1), matching the ε-dependence of matrix Hutch++.The result is stated for positive trace-class operators with a universal-constant high-probability bound.
- Truncated infHutch++: Basis-truncated infHutch++ converges in distribution to idealized infHutch++ for fixed sample budget, and its expectation converges to tr(A).The consistency argument also establishes convergence in probability under a coupling, but does not provide a finite-n Hutch++ error bound.
- Adaptive choice of sampling dimension: Adaptive basis selection stops when the newest block contributes little empirically, but this heuristic cannot certify the unresolved trace tail without additional decay assumptions.Several pilot samples are used to reduce stopping from a single atypically small realization.
4 Numerical examples
The numerical examples evaluate truncated infHutch and infHutch++ on integral operators, density-of-states estimation, and spectral filtering. Across these settings, the proposed truncations achieve comparable accuracy with lower-degree representations or smaller discretizations than ContHutch variants.
- Experimental setup: The experiments compare truncated infHutch and infHutch++ with ContHutch and ContHutch++ using operator applications as the primary cost metric.For Hutch++, the budget includes range-finding, residual estimation, and deterministic low-rank trace-term applications.
- 4.1 Integral kernel operators: Integral-kernel tests use Legendre-polynomial truncations and compare Helmholtz-like and sinc kernels against idealized infHutch benchmarks.The experiments also examine explicitly defined kernels and use idealized infHutch as a benchmark for truncated implementations.
- 4.1 Integral kernel operators: Selected truncated infHutch dimensions reach the idealized benchmark accuracy, with pilot-selected maximum degrees 20 for Helmholtz and 160 for sinc.The pilot procedure used m_pilot = 10 and basis tolerances interpolated between 10^-2 and 10^-3 across sample budgets.
- 4.1 Integral kernel operators: Comparable accuracy can require substantially lower Chebyshev degrees for truncated infHutch than for ContHutch, suggesting lower per-sample representation and operator-application cost.Increasing the Legendre degree beyond 100 did not further improve trace-estimation accuracy for either kernel.
- 4.2 Density-of-states: In density-of-states estimation, truncated infHutch++ reaches similar accuracy to ContHutch++ while ContHutch++ uses higher-degree samples, finer discretizations, and larger linear systems.Legendre samples of degree 512 reached the smallest attainable relative errors for each sample size, whereas degree 128 did not suffice.
- 4.3 Quadrature-based spectral filtering: For spectral filtering, the rational filter produces a positive trace-class approximation whose trace smoothly approximates the contour eigenvalue count.Eigenvalues well inside the contour receive weights near one, while those outside receive weights near zero; near-contour eigenvalues may be partially counted.
- 4.3 Quadrature-based spectral filtering: A second Dirac case has contour count 60 but full gap count 61, because an eigenvalue outside the contour is excluded by the filter.Counting that final eigenvalue would require moving the contour closer to the spectral edge and resolving a narrower transition region.
5 Conclusion
The paper establishes unbiased infinite-dimensional Girard–Hutchinson and Hutch++ estimators for positive trace-class operators, with practical truncations and favorable numerical behavior. Idealized infHutch++ achieves O(ε^-1) sample complexity, while truncated methods compare favorably with ContHutch++ in the reported experiments.
- Idealized estimators use target-operator-dependent Gaussian covariances, yielding unbiased infinite-dimensional versions of Girard–Hutchinson and Hutch++.
- O(ε^-1) sample complexity is proved for idealized infHutch++, matching the ε-dependence of matrix Hutch++.
- Truncated implementations sample from finite-dimensional subspaces; for infHutch++, the practical range finder restricts randomness while retaining operator images in the ambient space.
- In integral-operator experiments, truncated infHutch reaches idealized-benchmark accuracy with lower-degree Chebyshev representations than squared-exponential Gaussian-process samples.
- In density-of-states experiments, truncated infHutch++ achieves similar or better stochastic accuracy than ContHutch++ with substantially smaller internal chebfun discretizations.
6 Funding
The work was supported by Croatian Science Foundation grants and a project involving a Scientific Center of Excellence for Quantum and Complex Systems.
- The authors acknowledge support from Croatian Science Foundation grants IP-2025-02-3733 and IP-2022-10-5191.
- The acknowledged projects concern data-driven dynamical-systems identification, dimension reduction, optimal control, and model reduction.
- Additional support came from a project connected to the Scientific Center of Excellence for Quantum and Complex Systems and representations of Lie algebras.