Vecchia Gaussian Processes: Probabilistic Properties, Minimax Rates and Methodological Developments

Gaussian Processes (GPs) are widely used to model dependencies in spatial statistics and machine learning; however, exact inference is computationally intractable, with a time complexity of $O(n^3)$. Vecchia approximation allows scalable Bayesian inference of GPs by introducing sparsity in the spatial dependency structure characterized by a directed acyclic graph (DAG). Despite its practical popularity, the approach lacks rigorous theoretical foundations, and the choice of DAG structure remains an open question. In this paper, we systematically study Vecchia GPs as standalone stochastic processes, uncovering key probabilistic and statistical properties. For probabilistic properties, we prove that the conditional distributions of the Mat\'{e}rn GPs, as well as their Vecchia approximations, can be characterized by polynomial interpolations. This allows us to prove a series of results regarding small ball probabilities and Reproducing Kernel Hilbert Spaces (RKHSs) of Vecchia GPs. For the statistical methodology, we provide a principled guideline for selecting parent sets as norming sets with fixed cardinality, and we develop algorithms along these guidelines. In terms of theoretical guarantees, we establish posterior contraction rates for Vecchia GPs in the nonparametric regression model and show that minimax-optimal rates are attained under oracle rescaling or hierarchical Bayesian tuning. We demonstrate our theoretical results and methodology through numerical studies and provide an efficient implementation of our methods in C++, with interfaces in R.