
Abstract We introduce a Bayesian framework for geodesic regression on Grassmann manifolds. Since a geodesic on the manifold is fully determined by its initial position and velocity, we reformulate the regression from function estimation to inferring these two initial conditions. We perform this inference via a Grassmann-specific Lagrangian Hamiltonian Monte Carlo method that exploits the manifold’s Riemannian geometry. Solving the associated Lagrangian dynamics reduces to following the geodesic flow, yielding an efficient inference procedure. We analyze the asymptotic properties of the resulting posterior, providing theoretical guarantees, and confirm the computational efficiency and accuracy of the approach experimentally on complex Bayesian modelling tasks.
We introduce a mathematically new approach for quantization of vectorial signals. Namely, we are aimed at quantizing vectorial signals on input to a computational device that calculates some given function of the input. As opposed to the classical approach, which optimizes the quality of the quantized signal, without taking into account any further transformation of the latter by means of the computational device, here we optimize the quality of quantization on the output of this device. Moreover, we quantize components of the signal separately. This leads to a quantization problem qualitatively different from the classical one. We study existence of optimal quantizers (which is not at all obvious in this setting) and estimate the optimal cost for several classes of functions.
We study property testing for graphical models in a setting where an algorithm may adaptively query individual entries of the covariance matrix Sigma, and cost is measured by the number of queried entries. Under natural (strong) faithfulness assumptions, we design divide-and-conquer tests guided by balanced separators that operate on small submatrices and avoid global matrix inversion. Our first result is a tester for whether the underlying graph is a tree; with high probability it decides correctly using a subquadratic number of correlation queries, and for bounded maximum degree Delta the complexity is near-linear up to logarithmic factors. Our second result concerns graphs with a small separation number (hence small treewidth). We present two complementary procedures: a conditional descent test that never breaks when sn(G) <= k and, upon termination, returns an O(k log(n/k)) certificate; and a marginal descent test that, under a 'good run' condition, either certifies sn(G) <= 2k or proves sn(G) > k. The approach extends beyond the Gaussian case whenever reliable conditional-independence queries are available (e.g., non-paranormal models), yielding tests that access only a vanishing fraction of the entries of Sigma.
The main objective of this paper is to estimate optimally Sobol' indices at any order when a unique input/output i.i.d. sample is available. Our approach stands on three main ingredients: semi-parametric estimation theory, high-order kernel estimation (inspired by the paper Doksum & Samarov (1995, Ann. Stat., 23, 1443-1473)) and mirror-type transformations as introduced in Bertin et al. (2020, Electron. J. Stat., 14, 2198-2237) and Pujol (2020, arXiv:2205.03199) . We propose two different estimators. We prove that these estimators are asymptotically normal and efficient. Furthermore, we illustrate their numerical properties on standard examples.
Motivated by the classical work on finite noisy automata and by the recent work on broadcasting on grids, we introduce Gaussian variants of these models. These models are defined on graded posets. At time 0, all nodes begin with X-0. At time k >= 1, each node on layer k computes a combination of its inputs at layer k-1 with independent Gaussian noise added. When is it possible to recover X0 with non-vanishing correlation? We consider different notions of recovery including recovery from a single node, recovery from a bounded window and recovery from an unbounded window. Our main interest is in the following two models defined on grids. In the infinite model, layer k is the vertices of Z(d+1) whose sum of entries is k and for a vertex vat layer k >= 1, X-v = alpha (Xu + W-u,W-v), summed over all u on layer k-1 that differ from v exactly in one coordinate, and Wu,vare i.i.d. N (0, 1). We show that when alpha < 1/(d + 1), the correlation between Xv and X0 decays exponentially, and when alpha > 1/(d+ 1), the correlation is bounded away from 0. The critical case when alpha = 1/(d + 1) exhibits a phase transition in dimension, where X- v has non-vanishing correlation with X0 if and only if d >= 3. The same results hold for any bounded window. In the finite model, layer k is the vertices of Z(d+1) with non-negative entries with sum k. We identify the sub-critical and the super-critical regimes. In the sub-critical regime, the correlation decays to 0 for unbounded windows. In the super-critical regime, there exists for every t a convex combination of Xu on layer t whose correlation is bounded away from 0. We find that for the critical parameters, the correlation is vanishing in all dimensions and for unbounded window sizes.
This paper investigates the fundamental stability properties of the phase retrieval problem, which are critical for ensuring robust signal reconstruction across a wide range of applications. We present a rigorous analysis of the condition number $\beta _{\varPsi _{\boldsymbol A}}<^>{\ell _{p}}$ associated with the nonlinear mapping $\varPsi _{{\boldsymbol A}}(\boldsymbol{x}) = \bigl (\lvert \langle{\boldsymbol a}_{j}, {\boldsymbol x} angle vert <^>{2} \bigr )_{1 \le j \le m}$, evaluated with respect to the $\ell _{p}$ norm. Here, ${\boldsymbol a}_{j}\in{\mathbb R}<^>{n}$ or ${\mathbb C}<^>{n}$ denote the sensing vectors. A smaller condition number corresponds to enhanced stability against noise and perturbations. Our main contribution is the derivation of the first explicit and universal lower bounds on this condition number, which hold for all sensing vectors and thereby characterize intrinsic limitations on stability. Specifically, we establish that for the $\ell _{2}$ norm, the condition number $\beta _{\varPsi _{\boldsymbol A}}<^>{\ell _{2}}$ admits the following lower bounds: $$\begin{align*} \beta<^>{\ell_{2}}_{\varPsi_{\boldsymbol A}} \geq \begin{cases} \sqrt{3}, & ext{in the real case}, \\ 2, & ext{in the complex case}. \end{cases} \end{align*}$$Similarly, with respect to the $\ell _{1}$ norm, the condition number $\beta _{\varPsi _{\boldsymbol A}}<^>{\ell _{1}}$ asymptotically satisfies $$\begin{align*} \beta<^>{\ell_{1}}_{\varPsi_{\boldsymbol A}} \geq \begin{cases} \pi/2, & ext{in the real case}, \\ 2, & ext{in the complex case}. \end{cases} \end{align*}$$Importantly, we demonstrate that these fundamental bounds are asymptotically tight or exactly attained in significant instances. For example, in the real setting with $m imes 2$ sensing matrices, the harmonic frame achieves the optimal lower bound of $\sqrt{3}$ in the $\ell _{2}$ norm for all $m\geq 3$, thereby establishing itself as a provably optimal sensing matrix. By delineating these fundamental stability limits and identifying optimal measurement structures, our work provides profound insights into the intrinsic nature of phase retrieval and establishes a theoretical foundation for the development of more robust and efficient signal recovery algorithms.
We introduce a new methodology 'charcoal' for estimating the location of sparse changes in high-dimensional linear regression coefficients, without assuming that those coefficients are individually sparse. The procedure works by constructing different sketches (projections) of the design matrix at each time point so as to eliminate the possible dense nuisance parameters. The sequence of sketched design matrices is then compared against a single sketched response vector to form a sequence of test statistics whose behavior shows a surprising link to the well-known CUSUM statistics of univariate changepoint analysis. The procedure is computationally attractive, and strong theoretical guarantees are derived for its estimation accuracy. Simulations confirm that our methods perform well in extensive settings, and a real-world application to a large single-cell RNA sequencing dataset showcases the practical relevance.
This paper establishes sharp concentration inequalities for simple random tensors. Our theory unveils a phenomenon that arises only for asymmetric tensors of order $p \ge 3:$ when the effective ranks of the covariances of the component random variables lie on both sides of a critical threshold, an additional logarithmic factor emerges that is not present in sharp bounds for symmetric tensors. To establish our results, we develop empirical process theory for products of $p$ different function classes evaluated at $p$ different random variables, extending generic chaining techniques for quadratic and product empirical processes to higher-order settings.
We aim to approximate a continuously differentiable function $u:\mathbb{R}<^>{d} \rightarrow \mathbb{R}$ by a composition of functions $f\circ g$ where $g:\mathbb{R}<^>{d} \rightarrow \mathbb{R}<^>{m}$, $m\leq d$ and $f: \mathbb{R}<^>{m} \rightarrow \mathbb{R}$ are built in a two stage procedure. For a fixed $g$, we build $f$ using classical regression methods, involving evaluations of $u$. Recent works proposed to build a nonlinear $g$ by minimizing a loss function $\mathcal{J}(g)$ derived from Poincar & eacute; inequalities on manifolds, involving evaluations of the gradient of $u$. A problem is that minimizing $\mathcal{J}$ may be a challenging task. Hence in this work, we introduce new convex surrogates to $\mathcal{J}$. Leveraging concentration inequalities, we provide suboptimality results for a class of functions $g$, including polynomials, and a wide class of input probability measures. We investigate performances on different benchmarks for various training sample sizes. We show that our approach outperforms standard iterative methods for minimizing the training Poincar & eacute; inequality-based loss, often resulting in better approximation errors, especially for small training sets and $m=1$.
Randomized quasi-Monte Carlo (RQMC) methods estimate the mean of a random variable by sampling an integrand at $n$ equidistributed points. For scrambled digital nets, the resulting variance is typically $\tilde O(n<^>{-\theta })$, where $\theta \in [1,3]$ depends on the smoothness of the integrand and $\tilde O$ neglects logarithmic factors. While RQMC can be far more accurate than plain Monte Carlo (MC), it remains difficult to get confidence intervals on RQMC estimates. We investigate some empirical Bernstein confidence intervals (EBCIs) and hedged betting confidence intervals (HBCIs), both from Waudby-Smith and Ramdas (2024, J. Roy. Statist. Soc. B, 86, 1-27), when the random variable of interest is subject to known bounds. When there are $N$ integrand evaluations partitioned into $R$ independent replicates of $n=N/R$ RQMC points, and the RQMC variance is $\varTheta (n<^>{-\theta })$, then an oracle minimizing the width of a Bennett confidence interval would choose $n =\varTheta (N<^>{1/(\theta +1)})$. The resulting intervals have a width $\varTheta (N<^>{-\theta /(\theta +1)})$. Our empirical investigations had optimal values of $n$ grow slowly with $N$, HBCI intervals that were usually narrower than the EBCI ones and optimal values of $n$ for HBCI that were equal to or smaller than the ones for the oracle.
Given the ubiquity of streaming data, online algorithms have been widely used for parameter estimation, with second-order methods particularly standing out for their efficiency and robustness. In this paper, we study an online sketched Newton method that leverages a randomized sketching technique to perform an approximate Newton step in each iteration, thereby eliminating the computational bottleneck of second-order methods. While existing studies have established the asymptotic normality of sketched Newton methods, a consistent estimator of the limiting covariance matrix remains an open problem. We propose a fully online covariance matrix estimator that is constructed entirely from the Newton iterates and requires no matrix factorization. Compared to covariance estimators for first-order online methods, our estimator for second-order methods is batch-free. We establish the consistency and convergence rate of our estimator, and coupled with asymptotic normality results, we can then perform online statistical inference for the model parameters based on sketched Newton methods. We also discuss the extension of our estimator to constrained problems, and demonstrate its superior performance on regression problems as well as benchmark problems in the CUTEst set.
This paper demonstrates that when a shallow neural network with a Lipschitz continuous activation function is trained using either empirical or population risk to approximate a target function i.e. $r$ times continuously differentiable on $[0,1]<^>{d}$, the population risk may not decay at a rate faster than $t<^>{-\frac{4r}{d-2r}}$, where $t$ denotes the time parameter of the gradient flow dynamics. This result highlights the presence of the curse of dimensionality in the optimization computation required to achieve a desired accuracy. Instead of analyzing parameter evolution directly, the training dynamics are examined through the evolution of the parameter distribution under the 2-Wasserstein gradient flow. Furthermore, it is established that the curse of dimensionality persists when a locally Lipschitz continuous activation function is employed, where the Lipschitz constant in $[-x,x]$ is bounded by $O(x<^>\delta )$ for any $x \in \mathbb{R}$. In this scenario, the population risk is shown to decay at a rate no faster than $t<^>{-\frac{(4+2\delta )r}{d-2r}}$. Understanding how function smoothness influences the curse of dimensionality in neural network optimization theory is an important and underexplored direction that this work aims to address.
Score-based generative models, which transform noise into data by learning to reverse a diffusion process, have become a cornerstone of modern generative AI. This paper contributes to establishing theoretical guarantees for the probability flow ODE, a widely used diffusion-based sampler known for its practical efficiency. While a number of prior works address its general convergence theory, it remains unclear whether the probability flow ODE sampler can adapt to the low-dimensional structures commonly present in natural image data. We demonstrate that, with accurate score function estimation, the probability flow ODE sampler achieves a convergence rate of O(k/T) in total variation distance (ignoring logarithmic factors), where k is the intrinsic dimension of the target distribution and T is the number of iterations. This dimension-free convergence rate improves upon existing results that scale with the typically much larger ambient dimension, highlighting the ability of the probability flow ODE sampler to exploit intrinsic low-dimensional structures in the target distribution for faster sampling.
In this short note, we consider models of optimal Bayesian inference of finite-rank tensor products. We add to the model a linear channel parametrized by h. We show that at every interior differentiable point h of the free energy (associated with the model), the overlap concentrates at the gradient of the free energy and the minimum mean-square error converges to a related limit. In other words, the model is replica-symmetric at every differentiable point. At any signal-to-noise ratio, such points h form a full-measure set (hence h=0 belongs to the closure of these points). For a sufficiently low signal-to-noise ratio, we show that every interior point is a differentiable point.
Phase-only compressed sensing (PO-CS) concerns the recovery of sparse signals from the phases of complex measurements. Recent results show that sparse signals in the standard sphere $\mathbb{S}<^>{n-1}$ can be exactly recovered from complex Gaussian phases by a linearization procedure, which recasts PO-CS as linear compressed sensing and then applies (quadratically constrained) basis pursuit to obtain $ extbf{x}<^>\sharp$. This paper focuses on the instance optimality and robustness of $ extbf{x}<^>{\sharp }$. First, we strengthen the non-uniform instance optimality of Jacques and Feuillen (2021, IEEE Trans. Inf. Theory, 67, 4150-4161) to a uniform one over the entire signal space. We show the existence of some universal constant $C$ such that $\| extbf{x}<^>\sharp - extbf{x}\|_{2}\le Cs<^>{-1/2}\sigma _{\ell _{1}}( extbf{x},\varSigma <^>{n}_{s})$ holds for all $ extbf{x}$ in the unit Euclidean sphere, where $\sigma _{\ell _{1}}( extbf{x},\varSigma <^>{n}_{s})$ is the $\ell _{1}$ distance of $ extbf{x}$ to its closest $s$-sparse signal. This is achieved by showing that the new sensing matrices corresponding to all approximately sparse signals simultaneously satisfy restricted isometry property. Second, we investigate the estimator's robustness to noise and corruption. We show that dense noise with entries bounded by some small $ au _{0}$, appearing either prior or posterior to retaining the phases, increments $\| extbf{x}<^>\sharp - extbf{x}\|_{2}$ by $O( au _{0})$. This is near-optimal (up to log factors) for any algorithm. On the other hand, adversarial corruption, which changes an arbitrary $\zeta _{0}$-fraction of the measurements to any phase-only values, increments $\| extbf{x}<^>\sharp - extbf{x}\|_{2}$ by $O(\sqrt{\zeta _{0}\log (1/\zeta _{0})})$. We demonstrate the tightness of this result via a partial analysis under suboptimal noise parameter and numerical evidence, while showing that the impact of sparse corruption can be eliminated: to this end, we propose an extended linearization approach that can exactly recover $ extbf{x}$ from the corrupted phases. The developments are then combined to yield a robust instance optimal guarantee that resembles the standard one in linear compressed sensing.
We study the error introduced by entropy regularization in infinite-horizon discrete discounted Markov decision processes. We show that this error decreases exponentially in the inverse regularization strength, both in a weighted Kullback-Leibler divergence and in value with a problem-specific exponent. This is in contrast to previously known estimates, of the order $O(\tau )$, where $\tau$ is the regularization strength. We provide a lower bound that matches our upper bound up to a polynomial term, thereby characterizing the exponential convergence rate for entropy regularization. Our proof relies on the observation that the solutions of entropy-regularized Markov decision processes solve a gradient flow of the unregularized reward with respect to a Riemannian metric common in natural policy gradient methods. This correspondence allows us to identify the limit of this gradient flow as the generalized maximum entropy optimal policy, thereby characterizing the implicit bias of this gradient flow, which corresponds to a time-continuous version of the natural policy gradient method. We use our improved error estimates to show that for entropy-regularized natural policy gradient methods, the overall error decays exponentially in the square root of the number of iterations, improving over existing sublinear guarantees. Finally, we extend our analysis to settings beyond the entropy. In particular, we characterize the implicit bias regarding general convex potentials and their resulting generalized natural policy gradients.
This study proposes a novel method for estimation and hypothesis testing in high-dimensional single-index models. We address a common scenario where the sample size and the dimension of regression coefficients are large and comparable. Unlike previous approaches, which often overlook the estimation of the unknown link function, we introduce a new method for link function estimation. Leveraging the information from the estimated link function, we propose more efficient estimators that are better aligned with the underlying model. Furthermore, we rigorously establish the asymptotic normality of each coordinate of the estimator. This provides a valid construction of confidence intervals and p-values for any finite collection of coordinates. Numerical experiments validate our theoretical results.
Sparse signal recovery deals with finding the sparsest solution of an under-determined linear system $\boldsymbol{x} = \boldsymbol{Q}\boldsymbol{s}$. In this paper, we propose a novel greedy approach to addressing the challenges from such a problem. Such an approach is based on a characterization of solutions to the system, which allows us to work on the sparse recovery in the $\boldsymbol{s}$-space directly with a given measure. With $l_{2}$-based measure, an orthogonal matching pursuit (OMP)-type algorithm is proposed, which significantly outperforms the classical OMP algorithm in terms of recovery accuracy while maintaining comparable computational complexity. An $l_{1}$-based algorithm, denoted as $ ext{Alg}_{GL1}$, is derived. Such an algorithm significantly outperforms the classical basis pursuit algorithm. Combining with the compressive sampling match pursuit strategy for selecting atoms, a class of high-performance greedy algorithms is also derived. Extensive numerical simulations on both synthetic and image data are carried out, with which the superior performance of our proposed algorithms is demonstrated in terms of sparse recovery accuracy and robustness against numerical instability of the system matrix $\boldsymbol{Q}$ and disturbance in the measurement $\boldsymbol{x}$.
We consider $U$-statistics on row-column exchangeable matrices, which are arrays invariant to separate permutations of rows and columns and are common in bipartite data. Under the standard dissociation assumption, we develop a graph-indexed analogue of the Hoeffding decomposition tailored to row-column exchangeable dependence. We present a new decomposition based on orthogonal projections onto probability spaces generated by sets of Aldous-Hoover-Kallenberg variables. These sets are indexed by bipartite graphs, enabling the application of graph-theoretic concepts to describe the decomposition. This framework provides new insights into the characterization of $U$-statistics on row-column exchangeable matrices, particularly regarding their asymptotic behaviour, including in degenerate cases. Notably, the limit distribution depends only on specific terms in the decomposition, corresponding to non-zero components indexed by the smallest graphs, namely the principal support graphs. We show that the asymptotic behaviour of a $U$-statistic is characterized by the properties of its principal support graphs. The number of nodes in these graphs (the principal degree) dictates the convergence rate to the limit distribution, with degeneracy occurring if and only if this number is strictly greater than 1. Furthermore, when the principal support graphs are connected, the limit distribution is Gaussian, even in degenerate cases. Applications to network analysis illustrate these findings.
In this paper, we investigate the recovery of the sparse representation of data in general infinite-dimensional optimization problems regularized by convex functionals. We show that it is possible to define a suitable non-degeneracy condition on the minimal-norm dual certificate, extending the well-established non-degeneracy source condition (NDSC) associated with total variation regularized problems in the space of measures, as introduced in (Duval and Peyr\'e, FoCM, 15:1315-1355, 2015). In our general setting, we need to study how the dual certificate is acting, through the duality product, on the set of extreme points of the ball of the regularizer, seen as a metric space. This justifies the name Metric Non-Degenerate Source Condition (MNDSC). More precisely, we impose a second-order condition on the dual certificate, evaluated on curves with values in small neighbourhoods of a given collection of n extreme points. By assuming the validity of the MNDSC, together with the linear independence of the measurements on these extreme points, we establish that, for a suitable choice of regularization parameters and noise levels, the minimizer of the minimization problem is unique and is uniquely represented as a linear combination of n extreme points. The paper concludes by obtaining explicit formulations of the MNDSC for three problems of interest. First, we examine total variation regularized deconvolution problems, showing that the classical NDSC implies our MNDSC, and recovering a result similar to (Duval and Peyr\'e, FoCM, 15:1315-1355, 2015). Then, we consider 1-dimensional BV functions regularized with their BV-seminorm and pairs of measures regularized with their mutual 1-Wasserstein distance. In each case, we provide explicit versions of the MNDSC and formulate specific sparse representation recovery results.