Sublinear Time Low-Rank Approximation of Toeplitz Matrices

We present a sublinear time algorithm for computing a near optimal low-rank approximation to any positive semidefinite (PSD) Toeplitz matrix $T\in {\mathbb{R}}^{d\times d}$, given noisy access to its entries. In particular, given entrywise query access to $T+E$ for an arbitrary noise matrix $E\in {\mathbb{R}}^{d\times d}$, integer rank $k\le d$, and error parameter $\delta >0$, our algorithm runs in time $\text{poly}(k,\log (d\mathrm{/}\delta ))$ and outputs (in factored form) a Toeplitz matrix $\tilde{T}\in {\mathbb{R}}^{d\times d}$ with rank $\text{poly}(k,\log (d\mathrm{/}\delta ))$ satisfying, for some fixed constant $C$, \begin{equation*} |T-\widetilde{T}|_F \leq C \cdot \max{|E|_F,|T-T_k|_F} + δ\cdot |T|_F. \end{equation*} Here $\mathrm{\parallel }\cdot {\mathrm{\parallel }}_F$ is the Frobenius norm and $T_k$ is the best (not necessarily Toeplitz) rank-$k$ approximation to $T$ in the Frobenius norm, given by projecting $T$ onto its top $k$ eigenvectors. Our result has the following applications. When $E=0$, we obtain the first sublinear time near-relative-error low-rank approximation algorithm for PSD Toeplitz matrices, resolving the main open problem of Kapralov et al. SODA 23, whose algorithm had sublinear query complexity but exponential runtime. Our algorithm can also be applied to approximate the unknown Toeplitz covariance matrix of a multivariate Gaussian distribution, given sample access to this distribution, resolving an open question of Eldar et al. SODA 20. Our algorithm applies sparse Fourier transform techniques to recover a low-rank Toeplitz matrix using its Fourier structure. Our key technical contribution is the first polynomial time algorithm for \emph{discrete time off-grid} sparse Fourier recovery, which may be of independent interest.