We consider the task of fitting low-dimensional embeddings to high-dimensional data. In particular, we study the \(k\)-Euclidean Metric Violation problem (\(k-\textsf{EMV}\)), where the input is \(D \in \mathbb{R}_{\geqslant 0}^{\binom{n}{2}}\) and the goal is to find the closest vector \(X \in \mathbb{M}_k\), where \(\mathbb{M}_k \subset \mathbb{R}_{\geqslant 0}^{\binom{n}{2}}\) is the set of all \(k\)-dimensional Euclidean metrics on \(n\) points, and closeness is formulated as the following optimization problem, where \(\|\cdot\|\) is the entry-wise \(\ell_2\) norm: \(\mathsf{OPT}_{\textsf{EMV}} = \min_{X \in \mathbb{M}_k} \|D - X\|_2^2\). Cayton and Dasgupta [CD06] showed that this problem is NP-Hard, even when \(k = 1\). Dhamdhere [Dha04] obtained a \(O(\log(n))\)-approximation for \(1-\textsf{EMV}\) and leaves finding a PTAS for it as an open question (reiterated recently by Lee [Lee25]). Although \(k-\textsf{EMV}\) has been studied in the statistics community for over 70 years, under the name “multi-dimensional scaling,” there are no known efficient approximation algorithms for \(k \gt 1\), to the best of our knowledge.

A central challenge in quantum physics is to understand the structural properties of many-body systems, both in equilibrium and out of equilibrium. For classical systems, we have a unified perspective which connects structural properties of systems at thermal equilibrium to the Markov chain dynamics that mix to them. We lack such a perspective for quantum systems: there is no framework to translate the quantitative convergence of the Markovian evolution into strong structural consequences.

We give a polynomial-time algorithm for the problem of robustly estimating a mixture of k arbitrary Gaussians in ℝ d , for any fixed k , in the presence of a constant fraction of arbitrary corruptions. This resolves the main open problem in several previous works on algorithmic robust statistics, which addressed the special cases of robustly estimating (a) a single Gaussian, (b) a mixture of TV-distance separated Gaussians, and (c) a uniform mixture of two Gaussians. Our main tools are an efficient partial clustering algorithm that relies on the sum-of-squares method, and a novel tensor decomposition algorithm that allows errors in both Frobenius norm and low-rank terms.

Metric embeddings are a widely used method in algorithm design, where generally a "complex" metric is embedded into a simpler, lower-dimensional one. Historically, the theoretical computer science community has focused on bi-Lipschitz embeddings, which guarantee that every pairwise distance is approximately preserved. In contrast, alternative embedding objectives that are commonly used in practice avoid bi-Lipschitz distortion; yet these approaches have received comparatively less study in theory. In this paper, we focus on Multi-dimensional Scaling (MDS), where we are given a set of non-negative dissimilarities $\{d_{i,j}\}_{i,j\in[n]}$ over $n$ points, and the goal is to find an embedding $\{x_1,…,x_n\}\subset\mathbb{R}^k$ that minimizes \[ \mathrm{OPT} = \min_{x_1,…,x_n} \mathbb{E}_{i,j\in[n]} \left[ \left(1-\frac{\|x_i - x_j\|}{d_{i,j}}\right)^2 \right]. \]{Despite} its popularity, our theoretical understanding of MDS is extremely limited. Recently, Demaine et al. gave the first approximation algorithm with provable guarantees for this objective, which achieves an embedding in constant-dimensional Euclidean space with cost $\mathrm{OPT} + \epsilon$ in $n^2 \cdot 2^{\mathrm{poly}(\Delta/\epsilon)}$ time, where $\Delta$ is the aspect ratio of the input dissimilarities. For metrics that admit low-cost embeddings, $\Delta$ scales polynomially in $n$. In this work, we give the first approximation algorithm for MDS with quasi-polynomial dependency on $\Delta$: for constant-dimensional Euclidean space, we achieve a solution with cost $O(\log \Delta)\cdot \mathrm{OPT}^{\Omega(1)} + \epsilon$ in time $n^{O(1)} \cdot 2^{\mathrm{poly}\left(\frac{\log(\Delta)}{\epsilon}\right)}$. Our algorithms are based on a novel geometry-aware analysis of a conditional rounding of the Sherali-Adams LP hierarchy, allowing us to avoid the exponential dependency on the aspect ratio that would typically result from this rounding. \end{abstract}

STOC ’25, Prague, Czechia

We study the problem of finding a (pure) product state with optimal fidelity to an unknown $n$-qubit quantum state $ρ$, given copies of $ρ$. This is a basic instance of a fundamental question in quantum learning: is it possible to efficiently learn a simple approximation to an arbitrary state? We give an algorithm which finds a product state with fidelity $\varepsilon$-close to optimal, using $N = n^{\text{poly}(1/\varepsilon)}$ copies of $ρ$ and $\text{poly}(N)$ classical overhead. We further show that estimating the optimal fidelity is NP-hard for error $\varepsilon = 1/\text{poly}(n)$, showing that the error dependence cannot be significantly improved. For our algorithm, we build a carefully-defined cover over candidate product states, qubit by qubit, and then demonstrate that extending the cover can be reduced to approximate constrained polynomial optimization. For our proof of hardness, we give a formal reduction from polynomial optimization to finding the closest product state. Together, these results demonstrate a fundamental connection between these two seemingly unrelated questions. Building on our general approach, we also develop more efficient algorithms in three simpler settings: when the optimal fidelity exceeds $5/6$; when we restrict ourselves to a discrete class of product states; and when we are allowed to output a matrix product state.

A set of high dimensional points X$= \{x_{1},x_{2},\ldots,x_{n}\}\subseteq \mathbb{R}^{d}$in isotropic position is said to be$\delta$-anti concentrated if for every direction$v$, the fraction of points in$X$satisfying$\left|\left\langle x_i, v\right\rangle\right| \leqslant \delta$is at most$O(\delta)$. Motivated by applications to list-decodable learning and clustering, three recent works [7], [44], [71] considered the problem of constructing efficient certificates of anti-concentration in the average case, when the set of points X corresponds to samples from a Gaussian distribution. Their certificates played a crucial role in several subsequent works in algorithmic robust statistics on list-decodable learning and settling the robust learnability of arbitrary Gaussian mixtures. Unlike related efficient certificates of concentration properties that are known for wide class of distri-butions [52], the aforementioned approach has been limited only to rotationally invariant distributions (and their affine transformations) with the only prominent example being Gaussian distributions. This work presents a new (and arguably the most natural) formulation for anti- concentration. Using this formulation, we give quasi-polynomial time verifiable sum-of-squares certificates of anti-concentration that hold for a wide class of non-Gaussian distributions including anti-concentrated bounded product distributions and uniform distributions over$L_{p}$balls (and their affine transformations). Consequently, our method upgrades and extends results in algorithmic robust statistics e.g., list-decodable learning and clustering, to such distributions. As in the case of previous works, our certificates are also obtained via relaxations in the sum-of-squares hierarchy. However, the nature of our argument differs significantly from prior works that formulate anti-concentration as the non-negativity of an explicit polynomial. Our argument constructs a canonical integer program for anti-concentration and analysis a SoS relaxation of it, independent of the intended application. The explicit polynomials appearing in prior works can be seen as specific dual certificates to this program. From a technical standpoint, unlike existing works that explicitly construct sum-of-squares certificates, our argument relies on duality and analyzes a pseudo-expectation on large subsets of the input points that take a small value in some direction. Our analysis uses the method of polynomial reweightings to reduce the problem to analyzing only analytically dense or sparse directions.

We show that thermal states of local Hamiltonians are separable above a constant temperature. Specifically, for a local Hamiltonian$H$on a graph with degree$\mathfrak{g}$, its Gibbs state at inverse temperature$\beta$, denoted by$\rho=e^{-\beta H}/\text{tr}(e^{-\beta H})$, is a classical distribution over product states for all$\beta < 1/ (c \mathfrak{{g}})$, where$c$is a constant. This sudden death of thermal entanglement upends conventional wisdom about the presence of short-range quantum correlations in Gibbs states. Moreover, we show that we can efficiently sample from the distribution over product states. In particular, for any$\beta < 1/(c\mathfrak{g}^{3})$, we can prepare a state$\varepsilon$-close to$\rho$in trace distance with a depth-one quantum circuit and$\text{poly}(n)\log(1/\varepsilon)$classical overhead.11In independent and concurrent work, Rouzé, França, and Alhambra [37] obtain an efficient quantum algorithm for preparing high-temperature Gibbs states via a dissipative evolution.

We study the problem of Hamiltonian structure learning from real-time evolution: given the ability to apply$e^{-\mathrm{i}Ht}$for an unknown local Hamiltonian$H=\Sigma_{a=1}^{m}\lambda_{a}E_{a}$on$n$qubits, the goal is to recover$H$. This problem is already well-understood under the assumption that the interaction terms,$E_{a}$, are given, and only the interaction strengths,$\lambda_{a}$, are unknown. But how efficiently can we learn a local Hamiltonian without prior knowledge of its interaction structure? We present a new, general approach to Hamiltonian learning that not only solves the challenging structure learning variant, but also resolves other open questions in the area, all while achieving the gold standard of Heisenberg-limited scaling. In particular, our algorithm recovers the Hamiltonian to$\varepsilon$error with total evolution time$\mathcal{O}(\log(n)/\varepsilon)$, and has the following appealing properties: 1)It does not need to know the Hamiltonian terms; 2)It works beyond the short-range setting, extending to any Hamiltonian$H$where the sum of terms interacting with a qubit has bounded norm; 3)It evolves according to$H$in constant time$t$increments, thus achieving constant time resolution. As an application, we can also learn Hamiltonians exhibiting power-law decay up to accuracy$\varepsilon$with total evolution time beating the standard limit of$1/\varepsilon^{2}$.

The field of quantum machine learning (QML) produces many proposals for attaining quantum speedups for tasks in machine learning and data analysis. Such speedups can only manifest if classical algorithms for these tasks perform significantly slower than quantum ones. We study quantum-classical gaps in QML through the quantum singular value transformation (QSVT) framework. QSVT, introduced by Gilyén, Su, Low and Wiebe [GSLW19], unifies all major types of quantum speedup [MRTC21]; in particular, a wide variety of QML proposals are applications of QSVT on low-rank classical data. We challenge these proposals by providing a classical algorithm that matches the performance of QSVT in this regime up to a small polynomial overhead.

We study the problem of learning a local quantum Hamiltonian H given copies of its Gibbs state ρ = e−β H/(e−β H) at a known inverse temperature β>0. Anshu, Arunachalam, Kuwahara, and Soleimanifar gave an algorithm to learn a Hamiltonian on n qubits to precision with only polynomially many copies of the Gibbs state, but which takes exponential time. Obtaining a computationally efficient algorithm has been a major open problem, with prior work only resolving this in the limited cases of high temperature or commuting terms. We fully resolve this problem, giving a polynomial time algorithm for learning H to precision from polynomially many copies of the Gibbs state at any constant β > 0. Our main technical contribution is a new flat polynomial approximation to the exponential function, and a translation between multi-variate scalar polynomials and nested commutators. This enables us to formulate Hamiltonian learning as a polynomial system. We then show that solving a low-degree sum-of-squares relaxation of this polynomial system suffices to accurately learn the Hamiltonian.

We consider the problem of rank-1 low-rank approximation (LRA) in the matrix-vector product model under various Schatten norms:\begin{equation*}\min _{\|u\|_{2}=1}\left\|A\left(I-u u^{\top}\right)\right\|_{\mathcal{S}_{p}}\end{equation*}where $\|M\|_{\mathcal{S}_{p}}$ denotes the $\ell_{p}$ norm of the singular values of M. Given $\varepsilon\gt 0$, our goal is to output a unit vector v such that \begin{equation*}\left\|A\left(I-v v^{\top}\right)\right\|_{\mathcal{S}_{p}} \leqslant\left(1+\varepsilon\right) \min _{\|u\|_{2}=1}\left\|A\left(I-u u^{\top}\right)\right\|_{\mathcal{S}_{p}}\end{equation*}Our main result shows that Krylov methods (nearly) achieve the information-theoretically optimal1number of matrix-vector products for Spectral $(p=\infty)$, Frobenius $(p=2)$ and Nuclear $(p=1)$ LRA. In particular, for Spectral LRA, we show that any algorithm requires $\Omega\left(\log (n) / \varepsilon^{1 / 2}\right)$ matrix-vector products, exactly matching the upper bound obtained by Krylov methods [40]. Our lower bound addresses Open Question 1 in [59], providing evidence for the lack of progress on algorithms for Spectral LRA and resolves Open Question 1.2 in [5]. Next, we show that for any fixed constant p, i.e. $1 \leqslant p=O(1)$, there is an upper bound of $O\left(\log (1 / \varepsilon) / \varepsilon^{1 / 3}\right)$ matrix-vector products, implying that the complexity does not grow as a function of input size. This improves the $O\left(\log (n / \varepsilon) / \varepsilon^{1 / 3}\right)$ bound recently obtained in [5], and matches their $\Omega\left(1 / \varepsilon^{1 / 3}\right)$ lower bound, to a $\log (1 / \varepsilon)$ factor.1For Spectral LRA, the upper and lower bounds match up to a fixed universal constant. For Frobenius and Nuclear LRA, they match up to a $\log (1 / \varepsilon)$ factor.

Recently Chen and Poor initiated the study of learning mixtures of linear dynamical systems. While linear dynamical systems already have wide-ranging applications in modeling time-series data, using mixture models can lead to a better fit or even a richer understanding of underlying subpopulations represented in the data. In this work we give a new approach to learning mixtures of linear dynamical systems that is based on tensor decompositions. As a result, our algorithm succeeds without strong separation conditions on the components, and can be used to compete with the Bayes optimal clustering of the trajectories. Moreover our algorithm works in the challenging partially-observed setting. Our starting point is the simple but powerful observation that the classic Ho-Kalman algorithm is a relative of modern tensor decomposition methods for learning latent variable models. This gives us a playbook for how to extend it to work with more complicated generative models.

For any two point sets $A,B \subset \mathbb{R}^d$ of size up to $n$, the Chamfer distance from $A$ to $B$ is defined as $\texttt{CH}(A,B)=\sum_{a \in A} \min_{b \in B} d_X(a,b)$, where $d_X$ is the underlying distance measure (e.g., the Euclidean or Manhattan distance). The Chamfer distance is a popular measure of dissimilarity between point clouds, used in many machine learning, computer vision, and graphics applications, and admits a straightforward $O(d n^2)$-time brute force algorithm. Further, Chamfer distance is often used as a proxy for the more computationally demanding Earth-Mover (Optimal Transport) Distance. However, the \emph{quadratic} dependence on $n$ in the running time makes the naive approach intractable for large datasets. We overcome this bottleneck and present the first $(1+\epsilon)$-approximate algorithm for estimating Chamfer distance with a near-linear running time. Specifically, our algorithm runs in time $O(nd \log (n)/\epsilon^2)$ and is implementable. Our experiments demonstrate that it is both accurate and fast on large high-dimensional datasets. We believe that our algorithm will open new avenues for analyzing large high-dimensional point clouds. We also give evidence that if the goal is to report a $(1+\epsilon)$-approximate mapping from $A$ to $B$ (as opposed to just its value), then any sub-quadratic time algorithm is unlikely to exist.

Linear dynamical systems are the foundational statistical model upon which control theory is built. Both the celebrated Kalman filter and the linear quadratic regulator require knowledge of the system dynamics to provide analytic guarantees. Naturally, learning the dynamics of a linear dynamical system from linear measurements has been intensively studied since Rudolph Kalman's pioneering work in the 1960's. Towards these ends, we provide the first polynomial time algorithm for learning a linear dynamical system from a polynomial length trajectory up to polynomial error in the system parameters under essentially minimal assumptions; observability, controllability, and marginal stability. Our algorithm is built on a method of moments estimator to directly estimate Markov parameters from which the dynamics can be extracted. Furthermore we provide statistical lower bounds when our observability and controllability assumptions are violated.

We study iterative methods based on Krylov subspaces for low-rank approximation under any Schatten-p norm. Here, given access to a matrix A through matrix-vector products, an accuracy parameter є, and a target rank k, the goal is to find a rank-k matrix Z with orthonormal columns such that || A (I − Z Z⊤) || Sp ≤ (1+є)min U⊤U = Ik || A (I − U U⊤) || Sp, where || M ||Sp denotes the ℓp norm of the the singular values of M. For the special cases of p=2 (Frobenius norm) and p = ∞ (Spectral norm), Musco and Musco (NeurIPS 2015) obtained an algorithm based on Krylov methods that uses Õ(k/√є) matrix-vector products, improving on the naïve Õ(k/є) dependence obtainable by the power method, where Õ(·) suppresses poly(log(dk/є)) factors.

We give a polynomial-time algorithm for the problem of robustly estimating a mixture of k arbitrary Gaussians in ℝd, for any fixed k, in the presence of a constant fraction of arbitrary corruptions. This resolves the main open problem in several previous works on algorithmic robust statistics, which addressed the special cases of robustly estimating (a) a single Gaussian, (b) a mixture of TV-distance separated Gaussians, and (c) a uniform mixture of two Gaussians. Our main tools are an efficient partial clustering algorithm that relies on the sum-of-squares method, and a novel tensor decomposition algorithm that allows errors in both Frobenius norm and low-rank terms.

In list-decodable subspace recovery, the input is a collection of n points αn (for some α ≪ 1/2) of which are drawn i.i.d. from a distribution D with a isotropic rank r covariance Π∗ (the inliers) and the rest are arbitrary, potential adversarial outliers. The goal is to recover a O(1/α) size list of candidate covariances that contains a close to Π∗. Two recent independent works [56, 3] gave algorithms for this problem that work whenever D satisfies an algorithmic variant of anti-concentration condition (certifiable anticoncentration). The running time of both these algorithms, however, is and the error bounds on ‖Π – Π∗‖F grow with r (polynomially in r in [56] and logarithmically in [3]) that can be as large as Ω(d). In this work, we improve on these results on all three fronts: we obtain dimension-independent error in fixed-polynomial running time under less restrictive distributional assumptions. Specifically, we give a poly(1/α)dO(1) time algorithm that outputs a list containing a satisfying . Our result only needs certifiable hypercontractivity of degree 2 polynomials -a condition satisfied by a much broader family of distributions in contrast to certifiable anticoncentration. As a result, in addition to Gaussians, our algorithm applies to uniform distribution on the hypercube and q-ary cubes and arbitrary product distributions with subgaussian marginals. Prior work [56] had identified such distributions as potential hard examples as such distributions do not exhibit strong enough anti-concentration. When D satisfies certifiable anti-concentration, we obtain a stronger error guarantee of for any arbitrary η > 0 in dO(poly(1/α)+log(1/η)) time. The proof of the first result uses certifiable hypercontractivity of degree 2 polynomials to give a low-degree sum-of-squares proof of identifiability of the low dimensional structure in the presence of overwhelming fraction of outliers in input data. Our second result relies on a novel bootstrapping of the guarantees from the first with a new exponential error reduction mechanism within SoS along with certifiable anti-concentration. 1

We obtain robust and computationally efficient estimators for learning several linear models that achieve statistically optimal convergence rate under minimal distributional assumptions. Concretely, we assume our data is drawn from a k-hypercontractive distribution and an є-fraction is adversarially corrupted. We then describe an estimator that converges to the optimal least-squares minimizer for the true distribution at a rate proportional to є2−2/k, when the noise is independent of the covariates. We note that no such estimator was known prior to our work, even with access to unbounded computation. The rate we achieve is information-theoretically optimal and thus we resolve the main open question in Klivans, Kothari and Meka [COLT’18].

We study the Maximum Independent Set problem for geometric objects given in the data stream model. A set of geometric objects is said to be independent if the objects are pairwise disjoint. We consider geometric objects in one and two dimensions, i.e., intervals and disks. Let $α$ be the cardinality of the largest independent set. Our goal is to estimate $α$ in a small amount of space, given that the input is received as a one-pass stream. We also consider a generalization of this problem by assigning weights to each object and estimating $β$, the largest value of a weighted independent set. We initialize the study of this problem in the turnstile streaming model (insertions and deletions) and provide the first algorithms for estimating $α$ and $β$. For unit-length intervals, we obtain a $(2+ε)$-approximation to $α$ and $β$ in poly$(\frac{\log(n)}ε)$ space. We also show a matching lower bound. Combined with the $3/2$-approximation for insertion-only streams by Cabello and Perez-Lanterno [CP15], our result implies a separation between the insertion-only and turnstile model. For unit-radius disks, we obtain a $\left(\frac{8\sqrt{3}}π\right)$-approximation to $α$ and $β$ in poly$(\log(n), ε^{-1})$ space, which is closely related to the hexagonal circle packing constant. We provide algorithms for estimating $α$ for arbitrary-length intervals under a bounded intersection assumption and study the parameterized space complexity of estimating $α$ and $β$, where the parameter is the ratio of maximum to minimum interval length.

We study the problem of testing whether a matrix A ∈ \mathbbRn×nwith bounded entries ( ||A||∞ ≤ 1) is positive semidefinite (PSD), or ε-far in Euclidean distance from the PSD cone, meaning that , where B\succeq 0 denotes that B is PSD. Our main algorithmic contribution is a non-adaptive tester which distinguishes between these cases using only ~O(1/ε4) queries to the entries of A.11Throughout the paper, ~O(·) hides log(1/ε) factors. If instead of the Eucledian norm we considered the distance in spectral norm, we obtain the “ l∞-gap problem”, where A is either PSD or satisfies . For this related problem, we give a ~O(1/ε2) query tester, which we show is optimal up to log(1/ε) factors. Both our testers randomly sample a collection of principal sub-matrices and check whether these sub-matrices are PSD. Consequentially, our algorithms achieve one-sided error: whenever they output that A is not PSD, they return a certificate that A has negative eigenvalues. We complement our upper bound for PSD testing with Eucledian norm distance by giving a ~Ω(1/ε2) lower bound for any non-adaptive algorithm. Our lower bound construction is general, and can be used to derive lower bounds for a number of spectral testing problems. As an example of the applicability of our construction, we obtain a new ~Ω(1/ε4) sampling lower bound for testing the Schatten-1 norm with a εn1.5gap, extending a result of Balcan, Li, Woodruff, and Zhang [11]. In addition, our hard instance results in new sampling lower bounds for estimating the Ky-Fan Norm, and the cost of rank- k approximations, i.e. .

Recently, Musco and Woodruff (FOCS, 2017) showed that given an n×n positive semidefinite (PSD) matrix A, it is possible to compute a (1+ε-approximate relative-error low-rank approximation to A by querying Õ(nk/ε2.5) entries of A in time Õ(nk/ε2.5+nkω-1/ε2(ω-1)). They also showed that any relative-error low-rank approximation algorithm must query Ω(nk/ε) entries of A, this gap has since remained open. Our main result is to resolve this question by obtaining an optimal algorithm that queries Õ(nk/ε) entries of A and outputs a relative-error low-rank approximation in Õ(n·(k/ε)ω-1) time. Note, our running time improves that of Musco and Woodruff, and matches the information-theoretic lower bound if the matrix-multiplication exponent ω is 2. We then extend our techniques to negative-type distance matrices. Here, our input is a pair-wise distance matrix A corresponding to a point set P={x1, x2, ..., xn} such that Ai,j=||xi-xj||22. Bakshi and Woodruff (NeurIPS, 2018) showed a bi-criteria, relative-error low-rank approximation for negative-type metrics. Their algorithm queries Õ(nk/ε2.5) entries and outputs a rank-( k+4) matrix. We show that the bi-criteria guarantee is not necessary and obtain an Õ(nk/ε) query algorithm, which is optimal. Our algorithm applies to all distance matrices that arise from metrics satisfying negative-type inequalities, including l1,l2, spherical metrics, hypermetrics and effective resistances on a graph. We also obtain faster algorithms for ridge regression. Next, we introduce a new robust low-rank approximation model which captures PSD matrices that have been corrupted with noise. We assume that the Frobenius norm of the corruption is bounded. Here, we relax the notion of approximation to additive-error, since it is information-theoretically impossible to obtain a relative-error approximation in this setting. While a sample complexity lower bound precludes sublinear algorithms for arbitrary PSD matrices, we provide the first sublinear time and query algorithms when the corruption on the diagonal entries is bounded. As a special case, we show sample-optimal sublinear time algorithms for low-rank approximation of correlation matrices corrupted by noise.

We give the first outlier-robust efficient algorithm for clustering a mixture of k statistically separated d - dimensional Gaussians ( k-GMMs). Concretely, our algorithm takes input an ε-corrupted sample from a k-GMM and outputs an approximate clustering that misclassifies at most kO(k)(ε+η) fraction of the points whenever every pair of mixture components are separated by 1-exp(-poly(k/η)) in total variation distance. This is the statistically weakest possible notion of separation and allows, for e.g., clustering of mixtures with components with the same mean with covariances differing in a single unknown direction or separated in Frobenius distance. The running time of our algorithm is dpoly(k/η). Such results were not known prior to our work, even for k=2. More generally, our algorithms succeed for mixtures of any distribution that satisfies two well-studied analytic assumptions - sum-of-squares certifiable hypercontractivity and anti-concentration. As an immediate corollary, they extend to clustering mixtures of arbitrary affine transforms of the uniform distribution on the d-dimensional unit sphere. Even the information theoretic clusterability of separated distributions satisfying our analytic assumptions was not known and is likely to be of independent interest. Our algorithms build on the recent flurry of work relying on certifiable anti-concentration first introduced in [1], [2]. Our techniques expand the sum-of-squares toolkit to show robust certifiability of TV-separated Gaussian clusters in data. This involves giving a low-degree sum-of-squares proof of statements that relate parameter (i.e. mean and covariances) distance to total variation distance by relying only on hypercontractivity and anti-concentration.

We consider the following fundamental problem in the study of neural networks: given input examples $x \in \mathbb{R}^d$ and their vector-valued labels, as defined by an underlying generative neural network, recover the weight matrices of this network. We consider two-layer networks, mapping $\mathbb{R}^d$ to $\mathbb{R}^m$, with a single hidden layer and $k$ non-linear activation units $f(\cdot)$, where $f(x) = \max \{x , 0\}$ is the ReLU activation function. Such a network is specified by two weight matrices, $\mathbf{U}^* \in \mathbb{R}^{m \times k}, \mathbf{V}^* \in \mathbb{R}^{k \times d}$, such that the label of an example $x \in \mathbb{R}^{d}$ is given by $\mathbf{U}^* f(\mathbf{V}^* x)$, where $f(\cdot)$ is applied coordinate-wise. Given $n$ samples $x^1,…,x^n \in \mathbb{R}^d$ as a matrix $\mathbf{X} \in \mathbb{R}^{d \times n}$ and the label $\mathbf{U}^* f(\mathbf{V}^* \mathbf{X})$ of the network on these samples, our goal is to recover the weight matrices $\mathbf{U}^*$ and $\mathbf{V}^*$. More generally, our labels $\mathbf{U}^* f(\mathbf{V}^* \mathbf{X})$ may be corrupted by noise, and instead we observe $\mathbf{U}^* f(\mathbf{V}^* \mathbf{X}) + \mathbf{E}$ where $\mathbf{E}$ is some noise matrix. Even in this case, we may still be interested in recovering good approximations to the weight matrices $\mathbf{U}^*$ and $\mathbf{V}^*$. In this work, we develop algorithms and hardness results under varying assumptions on the input and noise. Although the problem is NP-hard even for $k=2$, by assuming Gaussian marginals over the input $\mathbf{X}$ we are able to develop polynomial time algorithms for the approximate recovery of $\mathbf{U}^*$ and $\mathbf{V}^*$. Perhaps surprisingly, in the noiseless case our algorithms recover $\mathbf{U}^*,\mathbf{V}^*$ \textit{exactly}, i.e. with no error, in \textit{strongly} polynomial time. To the best of the our knowledge, this is the first algorithm to accomplish exact recovery for the ReLU activation function. For the noisy case, we give the first polynomial time algorithm that approximately recovers the weights in the presence of mean-zero noise $\mathbf{E}$. Our algorithms generalize to a larger class of \textit{rectified} activation functions, $f(x) = 0$ when $x\leq 0$, and $f(x) > 0$ otherwise. Although our polynomial time results require $\mathbf{U}^*$ to have full column rank, we also give a fixed-parameter tractable algorithm (in $k$) when $\mathbf{U}^*$ does not have this property. Lastly, we give a fixed-parameter tractable algorithm for more arbitrary noise matrices $\mathbf{E}$, so long as they are independent of $\mathbf{X}$.

In this work, we study the k-median and k-means clustering problems when the data is distributed across many servers and can contain outliers. While there has been a lot of work on these problems for worst-case instances, we focus on gaining a finer understanding through the lens of beyond worst-case analysis. Our main motivation is the following: for many applications such as clustering proteins by function or clustering communities in a social network, there is some unknown target clustering, and the hope is that running a k-median or k-means algorithm will produce clusterings which are close to matching the target clustering. Worst-case results can guarantee constant factor approximations to the optimal k-median or k-means objective value, but not closeness to the target clustering. Our first result is a distributed algorithm which returns a near-optimal clustering assuming a natural notion of stability, namely, approximation stability [Awasthi and Balcan, 2014], even when a constant fraction of the data are outliers. The communication complexity is O~(sk+z) where s is the number of machines, k is the number of clusters, and z is the number of outliers. Next, we show this amount of communication cannot be improved even in the setting when the input satisfies various non-worst-case assumptions. We give a matching Omega(sk+z) lower bound on the communication required both for approximating the optimal k-means or k-median cost up to any constant, and for returning a clustering that is close to the target clustering in Hamming distance. These lower bounds hold even when the data satisfies approximation stability or other common notions of stability, and the cluster sizes are balanced. Therefore, Omega(sk+z) is a communication bottleneck, even for real-world instances.

Let $\PP=\{ p_1, p_2, \ldots p_n \}$ and $\QQ = \{ q_1, q_2 \ldots q_m \}$ be two point sets in an arbitrary metric space. Let $\AA$ represent the $m\times n$ pairwise distance matrix with $\AA_{i,j} = d(p_i, q_j)$. Such distance matrices are commonly computed in software packages and have applications to learning image manifolds, handwriting recognition, and multi-dimensional unfolding, among other things. In an attempt to reduce their description size, we study low rank approximation of such matrices. Our main result is to show that for any underlying distance metric $d$, it is possible to achieve an additive error low rank approximation in sublinear time. We note that it is provably impossible to achieve such a guarantee in sublinear time for arbitrary matrices $\AA$, and our proof exploits special properties of distance matrices. We develop a recursive algorithm based on additive projection-cost preserving sampling. We then show that in general, relative error approximation in sublinear time is impossible for distance matrices, even if one allows for bicriteria solutions. Additionally, we show that if $\PP = \QQ$ and $d$ is the squared Euclidean distance, which is not a metric but rather the square of a metric, then a relative error bicriteria solution can be found in sublinear time. Finally, we empirically compare our algorithm with the SVD and input sparsity time algorithms. Our algorithm is several hundred times faster than the SVD, and about $8$-$20$ times faster than input sparsity methods on real-world and and synthetic datasets of size $10^8$. Accuracy-wise, our algorithm is only slightly worse than that of the SVD (optimal) and input-sparsity time algorithms.

Pattern recognition is a vast field which has seen significant advances over the years. As the datasets under consideration grow larger and more comprehensive, using efficient techniques to process them becomes increasingly important. We present a versatile technique for the purpose of feature selection and extraction - Class Dependent Features (CDFs). CDFs identify the features innate to a class and extract them accordingly. The features thus extracted are relevant to the entire class and not just to the individual data item. This paper focuses on using CDFs to improve the accuracy of classification and at the same time control computational expense by tackling the curse of dimensionality. In order to demonstrate the generality of this technique, it is applied to two problem statements which have very little in common with each other - handwritten digit recognition and text categorization. It is found that for both problem statements, the accuracy is comparable to state-of-the-art results and the speed of the operation is considerably greater. Results are presented for Reuters-21578 and Web-KB datasets relating to text categorization and the MNIST and USPS datasets for handwritten digit recognition.

This paper presents a versatile technique for the purpose of feature selection and extraction - Class Dependent Features (CDFs). We use CDFs to improve the accuracy of classification and at the same time control computational expense by tackling the curse of dimensionality. In order to demonstrate the generality of this technique, it is applied to handwritten digit recognition and text categorization.