Cameron Musco

NeurIPS Conference 2025 Conference Paper

Distance Adaptive Beam Search for Provably Accurate Graph-Based Nearest Neighbor Search

• Yousef Al-Jazzazi
• Haya Diwan
• Jinrui Gou
• Cameron Musco
• Christopher Musco
• Torsten Suel

Nearest neighbor search is central in machine learning, information retrieval, and databases. For high-dimensional datasets, graph-based methods such as HNSW, DiskANN, and NSG have become popular thanks to their empirical accuracy and efficiency. These methods construct a directed graph over the dataset and perform beam search on the graph to find nodes close to a given query. While significant work has focused on practical refinements and theoretical understanding of graph-based methods, many questions remain. We propose a new distance-based termination condition for beam search to replace the commonly used condition based on beam width. We prove that, as long as the search graph is navigable, our resulting Adaptive Beam Search method is guaranteed to approximately solve the nearest-neighbor problem, establishing a connection between navigability and the performance of graph-based search. We also provide extensive experiments on our new termination condition for both navigable graphs and approximately navigable graphs used in practice, such as HNSW and Vamana graphs. We find that Adaptive Beam Search outperforms standard beam search over a range of recall values, data sets, graph constructions, and target number of nearest neighbors. It thus provides a simple and practical way to improve the performance of popular methods.

SODA Conference 2025 Conference Paper

Improved Spectral Density Estimation via Explicit and Implicit Deflation

• Rajarshi Bhattacharjee
• Rajesh Jayaram
• Cameron Musco
• Christopher Musco
• Archan Ray

TIST Journal 2025 Journal Article

Local Edge Dynamics and Opinion Polarization

• Nikita Bhalla
• Adam Lechowicz
• Cameron Musco

The proliferation of social media platforms, recommender systems, and their joint societal impacts have prompted significant interest in opinion formation and evolution within social networks. We study how local edge dynamics can drive opinion polarization. In particular, we introduce a variant of the classic Friedkin-Johnsen opinion dynamics, augmented with a simple time-evolving network model. Edges are iteratively added or deleted according to simple rules, modeling decisions based on individual preferences and network recommendations. Via simulations on synthetic and real-world graphs, we find that the combined presence of two dynamics gives rise to high polarization: (1) confirmation bias —i.e., the preference for nodes to connect to other nodes with similar expressed opinions and (2) friend-of-friend link recommendations, which encourage new connections between closely connected nodes. We show that our model is tractable to theoretical analysis, which helps explain how these local dynamics erode connectivity across opinion groups, affecting polarization and a related measure of disagreement across edges. Finally, we validate our model against real-world data, showing that our edge dynamics drive the structure of arbitrary graphs, including random graphs, to more closely resemble real social networks. Our code and supplemental materials are available at https://github.com/adamlechowicz/opinion-polarization/. 1

ICML Conference 2025 Conference Paper

Near-Optimal Consistency-Robustness Trade-Offs for Learning-Augmented Online Knapsack Problems

• Mohammadreza Daneshvaramoli
• Helia Karisani
• Adam Lechowicz
• Bo Sun 0004
• Cameron Musco
• Mohammad Hajiesmaili

This paper introduces a family of learning-augmented algorithms for online knapsack problems that achieve near Pareto-optimal consistency-robustness trade-offs through a simple combination of trusted learning-augmented and worst-case algorithms. Our approach relies on succinct, practical predictions—single values or intervals estimating the minimum value of any item in an offline solution. Additionally, we propose a novel fractional-to-integral conversion procedure, offering new insights for online algorithm design.

SODA Conference 2025 Conference Paper

Near-optimal hierarchical matrix approximation from matrix-vector products

• Tyler Chen
• Feyza Duman Keles
• Diana Halikias
• Cameron Musco
• Christopher Musco
• David Persson

NeurIPS Conference 2024 Conference Paper

Efficient and Private Marginal Reconstruction with Local Non-Negativity

• Brett Mullins
• Miguel Fuentes
• Yingtai Xiao
• Daniel Kifer
• Cameron Musco
• Daniel Sheldon

Differential privacy is the dominant standard for formal and quantifiable privacy and has been used in major deployments that impact millions of people. Many differentially private algorithms for query release and synthetic data contain steps that reconstruct answers to queries from answers to other queries that have been measured privately. Reconstruction is an important subproblem for such mechanisms to economize the privacy budget, minimize error on reconstructed answers, and allow for scalability to high-dimensional datasets. In this paper, we introduce a principled and efficient postprocessing method ReM (Residuals-to-Marginals) for reconstructing answers to marginal queries. Our method builds on recent work on efficient mechanisms for marginal query release, based on making measurements using a residual query basis that admits efficient pseudoinversion, which is an important primitive used in reconstruction. An extension GReM-LNN (Gaussian Residuals-to-Marginals with Local Non-negativity) reconstructs marginals under Gaussian noise satisfying consistency and non-negativity, which often reduces error on reconstructed answers. We demonstrate the utility of ReM and GReM-LNN by applying them to improve existing private query answering mechanisms.

NeurIPS Conference 2024 Conference Paper

Gaussian Process Bandits for Top-k Recommendations

• Mohit Yadav
• Daniel Sheldon
• Cameron Musco

Algorithms that utilize bandit feedback to optimize top-k recommendations are vital for online marketplaces, search engines, and content platforms. However, the combinatorial nature of this problem poses a significant challenge, as the possible number of ordered top-k recommendations from $n$ items grows exponentially with $k$. As a result, previous work often relies on restrictive assumptions about the reward or bandit feedback models, such as assuming that the feedback discloses rewards for each recommended item rather than a single scalar feedback for the entire set of top-k recommendations. We introduce a novel contextual bandit algorithm for top-k recommendations, leveraging a Gaussian process with a Kendall kernel to model the reward function. Our algorithm requires only scalar feedback from the top-k recommendations and does not impose restrictive assumptions on the reward structure. Theoretical analysis confirms that the proposed algorithm achieves sub-linear regret in relation to the number of rounds and arms. Additionally, empirical results using a bandit simulator demonstrate that the proposed algorithm outperforms other baselines across various scenarios.

NeurIPS Conference 2024 Conference Paper

Navigable Graphs for High-Dimensional Nearest Neighbor Search: Constructions and Limits

• Haya Diwan
• Jinrui Gou
• Cameron Musco
• Christopher Musco
• Torsten Suel

There has been significant recent interest in graph-based nearest neighbor search methods, many of which are centered on the construction of (approximately) "navigable" graphs over high-dimensional point sets. A graph is navigable if we can successfully move from any starting node to any target node using a greedy routing strategy where we always move to the neighbor that is closest to the destination according to the given distance function. The complete graph is obviously navigable for any point set, but the important question for applications is if sparser graphs can be constructed. While this question is fairly well understood in low-dimensions, we establish some of the first upper and lower bounds for high-dimensional point sets. First, we give a simple and efficient way to construct a navigable graph with average degree $O(\sqrt{n \log n })$ for any set of $n$ points, in any dimension, for any distance function. We compliment this result with a nearly matching lower bound: even under the Euclidean metric in $O(\log n)$ dimensions, a random point set has no navigable graph with average degree $O(n^{\alpha})$ for any $\alpha < 1/2$. Our lower bound relies on sharp anti-concentration bounds for binomial random variables, which we use to show that the {near-neighborhoods} of a set of random points do not overlap significantly, forcing any navigable graph to have many edges.

NeurIPS Conference 2024 Conference Paper

Nearly Optimal Approximation of Matrix Functions by the Lanczos Method

• Noah Amsel
• Tyler Chen
• Anne Greenbaum
• Cameron Musco
• Christopher Musco

Approximating the action of a matrix function $f(\vec{A})$ on a vector $\vec{b}$ is an increasingly important primitive in machine learning, data science, and statistics, with applications such as sampling high dimensional Gaussians, Gaussian process regression and Bayesian inference, principle component analysis, and approximating Hessian spectral densities. Over the past decade, a number of algorithms enjoying strong theoretical guarantees have been proposed for this task. Many of the most successful belong to a family of algorithms called Krylov subspace methods. Remarkably, a classic Krylov subspace method, called the Lanczos method for matrix functions (Lanczos-FA), frequently outperforms newer methods in practice. Our main result is a theoretical justification for this finding: we show that, for a natural class of rational functions, Lanczos-FA matches the error of the best possible Krylov subspace method up to a multiplicative approximation factor. The approximation factor depends on the degree of $f(x)$'s denominator and the condition number of $\vec{A}$, but not on the number of iterations $k$. Our result provides a strong justification for the excellent performance of Lanczos-FA, especially on functions that are well approximated by rationals, such as the matrix square root.

ICML Conference 2024 Conference Paper

On the Role of Edge Dependency in Graph Generative Models

• Sudhanshu Chanpuriya
• Cameron Musco
• Konstantinos Sotiropoulos
• Charalampos E. Tsourakakis

We investigate the trade-off between the representation power of graph generative models and model overlap, i. e. , the degree to which the model generates diverse outputs versus regurgitating its training data. In particular, we delineate a nested hierarchy of graph generative models categorized into three levels of complexity: edge independent, node independent, and arbitrarily dependent models. This hierarchy encapsulates a wide range of prevalent methods. We derive theoretical bounds on the number of triangles and other short-length cycles producible by each level of the hierarchy, finding that more complex dependency structure allows an improved trade-off between representation power and overlap. We provide instances demonstrating the asymptotic optimality of our bounds. Furthermore, we introduce new generative models for each of the three hierarchical levels, leveraging dense subgraph discovery. Our evaluation, conducted on real-world datasets, focuses on assessing the output quality and overlap of our proposed models in comparison to other popular models. Our results indicate that our simple, interpretable models provide competitive baselines to popular generative models. Through this investigation, we offer a structured and robust evaluation scheme, thereby facilitating the development of models capable of generating accurate and edge-diverse graphs.

SODA Conference 2024 Conference Paper

On the Unreasonable Effectiveness of Single Vector Krylov Methods for Low-Rank Approximation

• Raphael A. Meyer
• Cameron Musco
• Christopher Musco

SODA Conference 2024 Conference Paper

Sublinear Time Low-Rank Approximation of Toeplitz Matrices

• Cameron Musco
• Kshiteej Sheth

ICLR Conference 2023 Conference Paper

Direct Embedding of Temporal Network Edges via Time-Decayed Line Graphs

• Sudhanshu Chanpuriya
• Ryan A. Rossi
• Sungchul Kim
• Tong Yu 0001
• Jane Hoffswell
• Nedim Lipka
• Shunan Guo
• Cameron Musco

Temporal networks model a variety of important phenomena involving timed interactions between entities. Existing methods for machine learning on temporal networks generally exhibit at least one of two limitations. First, many methods assume time to be discretized, so if the time data is continuous, the user must determine the discretization and discard precise time information. Second, edge representations can only be calculated indirectly from the nodes, which may be suboptimal for tasks like edge classification. We present a simple method that avoids both shortcomings: construct the line graph of the network, which includes a node for each interaction, and weigh the edges of this graph based on the difference in time between interactions. From this derived graph, edge representations for the original network can be computed with efficient classical methods. The simplicity of this approach facilitates explicit theoretical analysis: we can constructively show the effectiveness of our method's representations for a natural synthetic model of temporal networks. Empirical results on real-world networks demonstrate our method's efficacy and efficiency on both link classification and prediction.

NeurIPS Conference 2023 Conference Paper

Exact Representation of Sparse Networks with Symmetric Nonnegative Embeddings

• Sudhanshu Chanpuriya
• Ryan Rossi
• Anup B. Rao
• Tung Mai
• Nedim Lipka
• Zhao Song
• Cameron Musco

Graph models based on factorization of the adjacency matrix often fail to capture network structures related to links between dissimilar nodes (heterophily). We introduce a novel graph factorization model that leverages two nonnegative vectors per node to interpretably account for links between both similar and dissimilar nodes. We prove that our model can exactly represent any graph with low arboricity, a property that many real-world networks satisfy; our proof also applies to related models but has much greater scope than the closest prior bound, which is based on low max degree. Our factorization also has compelling properties besides expressiveness: due to its symmetric structure and nonnegativity, fitting the model inherently finds node communities, and the model's link predictions can be interpreted in terms of these communities. In experiments on real-world networks, we demonstrate our factorization's effectiveness on a variety of tasks, including community detection and link prediction.

NeurIPS Conference 2023 Conference Paper

Finite Population Regression Adjustment and Non-asymptotic Guarantees for Treatment Effect Estimation

• Mehrdad Ghadiri
• David Arbour
• Tung Mai
• Cameron Musco
• Anup B. Rao

The design and analysis of randomized experiments is fundamental to many areas, from the physical and social sciences to industrial settings. Regression adjustment is a popular technique to reduce the variance of estimates obtained from experiments, by utilizing information contained in auxiliary covariates. While there is a large literature within the statistics community studying various approaches to regression adjustment and their asymptotic properties, little focus has been given to approaches in the finite population setting with non-asymptotic accuracy bounds. Further, prior work typically assumes that an entire population is exposed to an experiment, whereas practitioners often seek to minimize the number of subjects exposed to an experiment, for ethical and pragmatic reasons. In this work, we study the problems of estimating the sample mean, individual treatment effects, and average treatment effect with regression adjustment. We propose approaches that use techniques from randomized numerical linear algebra to sample a subset of the population on which to perform an experiment. We give non-asymptotic accuracy bounds for our methods and demonstrate that they compare favorably with prior approaches.

SODA Conference 2023 Conference Paper

Near-Linear Sample Complexity for L p Polynomial Regression

• Raphael A. Meyer
• Cameron Musco
• Christopher Musco
• David P. Woodruff
• Samson Zhou

We study L p polynomial regression. Given query access to a function f: [−1, 1]→ℝ, the goal is to find a degree d polynomial q̂ such that, for a given parameter ε > 0 Here || · || p is the L p norm, ‖ g ‖ p = (∫ 1 −1 | g(t) | p dt ) 1/ p. We show that querying f at points randomly drawn from the Chebyshev measure on [-1, 1] is a near-optimal strategy for polynomial regression in all L p norms. In particular, to find q̂, it suffices to sample points from [-1, 1] with probabilities proportional to this measure. While the optimal sample complexity for polynomial regression was well understood for L 2 and L ∞, our result is the first that achieves sample complexity linear in d and error (1 + ε) for other values of p without any assumptions. Our result requires two main technical contributions. The first concerns p ≤ 2, for which we provide explicit bounds on the L p Lewis weight function of the infinite linear operator underlying polynomial regression. Using tools from the orthogonal polynomial literature, we show that this function is bounded by the Chebyshev density. Our second key contribution is to take advantage of the structure of polynomials to reduce the p > 2 case to the p ≤ 2 case. By doing so, we obtain a better sample complexity than what is possible for general p -norm linear regression problems, for which Ω (d p/2 ) samples are required.

NeurIPS Conference 2023 Conference Paper

No-regret Algorithms for Fair Resource Allocation

• Abhishek Sinha
• Ativ Joshi
• Rajarshi Bhattacharjee
• Cameron Musco
• Mohammad Hajiesmaili

We consider a fair resource allocation problem in the no-regret setting against an unrestricted adversary. The objective is to allocate resources equitably among several agents in an online fashion so that the difference of the aggregate $\alpha$-fair utilities of the agents achieved by an optimal static clairvoyant allocation and the online policy grows sublinearly with time. The problem inherits its difficulty from the non-separable nature of the global $\alpha$-fairness function. Previously, it was shown that no online policy could achieve a sublinear standard regret in this problem. In this paper, we propose an efficient online resource allocation policy, called Online Fair Allocation ($\texttt{OFA}$), that achieves sublinear $c_\alpha$-approximate regret with approximation factor $c_\alpha=(1-\alpha)^{-(1-\alpha)}\leq 1. 445, $ for $0\leq \alpha < 1$. Our upper bound on the $c_\alpha$-regret for this problem exhibits a surprising \emph{phase transition} phenomenon -- transitioning from a power-law to a constant at the critical exponent $\alpha=\frac{1}{2}. $ Our result also resolves an open problem in designing an efficient no-regret policy for the online job scheduling problem in certain parameter regimes. Along the way, we introduce new algorithmic and analytical techniques, including greedy estimation of the future gradients for non-additive global reward functions and bootstrapping second-order regret bounds, which may be of independent interest.

SODA Conference 2023 Conference Paper

Toeplitz Low-Rank Approximation with Sublinear Query Complexity

• Michael Kapralov
• Hannah Lawrence
• Mikhail Makarov
• Cameron Musco
• Kshiteej Sheth

We present a sublinear query algorithm for outputting a near-optimal low-rank approximation to any positive semidefinite Toeplitz matrix T ∈ ℝ d×d. In particular, for any integer rank k ≤ d and ε, δ > 0, our algorithm makes Õ ( k 2 · log(1/δ) · poly(1/ε)) queries to the entries of T and outputs a rank Õ ( k · log(1/δ)/ε) matrix d×d such that || T – || F ≤ (1 + ε) · || T - T k || F + δ|| Τ || F. Here, || · || F is the Frobenius norm and T k is the optimal rank- k approximation to T, given by projection onto its top k eigenvectors. Õ (·) hides polylog( d ) factors. Our algorithm is structure-preserving, in that the approximation is also Toeplitz. A key technical contribution is a proof that any positive semidefinite Toeplitz matrix in fact has a near-optimal low-rank approximation which is itself Toeplitz. Surprisingly, this basic existence result was not previously known. Building on this result, along with the well-established off-grid Fourier structure of Toeplitz matrices [Cybenko'82], we show that Toeplitz with near optimal error can be recovered with a small number of random queries via a leverage-score-based off-grid sparse Fourier sampling scheme.

FOCS Conference 2022 Conference Paper

Active Linear Regression for ℓp Norms and Beyond

• Cameron Musco
• Christopher Musco
• David P. Woodruff
• Taisuke Yasuda 0002

We study active sampling algorithms for linear regression, which aim to query only a small number of entries of a target vector and output a near minimizer to the objective function. For $\ell_{p}$ norm regression for any $ 0\lt p \lt \infty$, we give an algorithm based on Lewis weight sampling which outputs $\mathrm{a}(1+\epsilon)$-approximate solution using just $\tilde{O}(d/\epsilon^{2})$ queries to b for $p\in(0, 1)$, $\tilde{O}$ $(d/\epsilon)$ queries for $p\in(1, 2)$, and $\tilde{O}$ $(d^{p/2}/\epsilon^{p})$ queries for $p\in(2, \ \infty)$. For $p\in(0, 2)$, our bounds are optimal up to logarithmic factors, thus settling the query complexity for this range of p. For $p\in(2, \ \infty)$, our dependence on d is optimal, while our dependence on $\epsilon$ is off by at most a single $\epsilon$ factor, up to logarithmic factors. Our result resolves an open question of Chen and Dereziński, who gave near optimal bounds for the $\ell_{1}$ norm, but required at least $d^{2}/\epsilon^{2}$ samples for $\ell_{p}$ regression with $p\in(1, 2)$, and gave no bounds for $p\in(2, \ \infty)$ or $p\in(0, 1)$. We also provide the first total sensitivity upper bound for loss functions with at most degree p polynomial growth. This improves a recent result of Tukan, Maalouf, and Feldman. By combining this with our techniques for $\ell_{p}$ regression, we obtain the first active regression algorithms for such loss functions, including the important cases of the Tukey and Huber losses. This answers another question of Chen and Dereziński. Our sensitivity bounds also give improvements to a variety of previous results using sensitivity sampling, including Orlicz norm subspace embeddings, robust subspace approximation, and dimension reduction for smoothed p-norms. Finally, our active sampling results give the first sublinear time algorithms for Kronecker product regression under every $\ell_{p}$ norm. Previous results required reading the entire b vector in the kernel feature space. 1 1 Extended abstract; full version available at https: //arxiv. org/abs/2111. 04888.

ICLR Conference 2022 Conference Paper

Fast Regression for Structured Inputs

• Raphael A. Meyer
• Cameron Musco
• Christopher Musco
• David P. Woodruff
• Samson Zhou

We study the $\ell_p$ regression problem, which requires finding $\mathbf{x}\in\mathbb R^{d}$ that minimizes $\|\mathbf{A}\mathbf{x}-\mathbf{b}\|_p$ for a matrix $\mathbf{A}\in\mathbb R^{n \times d}$ and response vector $\mathbf{b}\in\mathbb R^{n}$. There has been recent interest in developing subsampling methods for this problem that can outperform standard techniques when $n$ is very large. However, all known subsampling approaches have run time that depends exponentially on $p$, typically, $d^{\mathcal{O}(p)}$, which can be prohibitively expensive. We improve on this work by showing that for a large class of common \emph{structured matrices}, such as combinations of low-rank matrices, sparse matrices, and Vandermonde matrices, there are subsampling based methods for $\ell_p$ regression that depend polynomially on $p$. For example, we give an algorithm for $\ell_p$ regression on Vandermonde matrices that runs in time $\mathcal{O}(n\log^3 n+(dp^2)^{0.5+\omega}\cdot\text{polylog}\,n)$, where $\omega$ is the exponent of matrix multiplication. The polynomial dependence on $p$ crucially allows our algorithms to extend naturally to efficient algorithms for $\ell_\infty$ regression, via approximation of $\ell_\infty$ by $\ell_{\mathcal{O}(\log n)}$. Of practical interest, we also develop a new subsampling algorithm for $\ell_p$ regression for arbitrary matrices, which is simpler than previous approaches for $p \ge 4$.

NeurIPS Conference 2022 Conference Paper

Kernel Interpolation with Sparse Grids

• Mohit Yadav
• Daniel R. Sheldon
• Cameron Musco

Structured kernel interpolation (SKI) accelerates Gaussian processes (GP) inference by interpolating the kernel covariance function using a dense grid of inducing points, whose corresponding kernel matrix is highly structured and thus amenable to fast linear algebra. Unfortunately, SKI scales poorly in the dimension of the input points, since the dense grid size grows exponentially with the dimension. To mitigate this issue, we propose the use of sparse grids within the SKI framework. These grids enable accurate interpolation, but with a number of points growing more slowly with dimension. We contribute a novel nearly linear time matrix-vector multiplication algorithm for the sparse grid kernel matrix. We also describe how sparse grids can be combined with an efficient interpolation scheme based on simplicial complexes. With these modifications, we demonstrate that SKI can be scaled to higher dimensions while maintaining accuracy, for both synthetic and real datasets.

NeurIPS Conference 2022 Conference Paper

Modeling Transitivity and Cyclicity in Directed Graphs via Binary Code Box Embeddings

• Dongxu Zhang
• Michael Boratko
• Cameron Musco
• Andrew McCallum

Modeling directed graphs with differentiable representations is a fundamental requirement for performing machine learning on graph-structured data. Geometric embedding models (e. g. hyperbolic, cone, and box embeddings) excel at this task, exhibiting useful inductive biases for directed graphs. However, modeling directed graphs that both contain cycles and some element of transitivity, two properties common in real-world settings, is challenging. Box embeddings, which can be thought of as representing the graph as an intersection over some learned super-graphs, have a natural inductive bias toward modeling transitivity, but (as we prove) cannot model cycles. To this end, we propose binary code box embeddings, where a learned binary code selects a subset of graphs for intersection. We explore several variants, including global binary codes (amounting to a union over intersections) and per-vertex binary codes (allowing greater flexibility) as well as methods of regularization. Theoretical and empirical results show that the proposed models not only preserve a useful inductive bias of transitivity but also have sufficient representational capacity to model arbitrary graphs, including graphs with cycles.

NeurIPS Conference 2022 Conference Paper

Sample Constrained Treatment Effect Estimation

• Raghavendra Addanki
• David Arbour
• Tung Mai
• Cameron Musco
• Anup Rao

Treatment effect estimation is a fundamental problem in causal inference. We focus on designing efficient randomized controlled trials, to accurately estimate the effect of some treatment on a population of $n$ individuals. In particular, we study \textit{sample-constrained treatment effect estimation}, where we must select a subset of $s \ll n$ individuals from the population to experiment on. This subset must be further partitioned into treatment and control groups. Algorithms for partitioning the entire population into treatment and control groups, or for choosing a single representative subset, have been well-studied. The key challenge in our setting is jointly choosing a representative subset and a partition for that set. We focus on both individual and average treatment effect estimation, under a linear effects model. We give provably efficient experimental designs and corresponding estimators, by identifying connections to discrepancy minimization and leverage-score-based sampling used in randomized numerical linear algebra. Our theoretical results obtain a smooth transition to known guarantees when $s$ equals the population size. We also empirically demonstrate the performance of our algorithms.

NeurIPS Conference 2022 Conference Paper

Simplified Graph Convolution with Heterophily

• Sudhanshu Chanpuriya
• Cameron Musco

Recent work has shown that a simple, fast method called Simple Graph Convolution (SGC) (Wu et al. , 2019), which eschews deep learning, is competitive with deep methods like graph convolutional networks (GCNs) (Kipf & Welling, 2017) in common graph machine learning benchmarks. The use of graph data in SGC implicitly assumes the common but not universal graph characteristic of homophily, wherein nodes link to nodes which are similar. Here we confirm that SGC is indeed ineffective for heterophilous (i. e. , non-homophilous) graphs via experiments on synthetic and real-world datasets. We propose Adaptive Simple Graph Convolution (ASGC), which we show can adapt to both homophilous and heterophilous graph structure. Like SGC, ASGC is not a deep model, and hence is fast, scalable, and interpretable; further, we can prove performance guarantees on natural synthetic data models. Empirically, ASGC is often competitive with recent deep models at node classification on a benchmark of real-world datasets. The SGC paper questioned whether the complexity of graph neural networks is warranted for common graph problems involving homophilous networks; our results similarly suggest that, while deep learning often achieves the highest performance, heterophilous structure alone does not necessitate these more involved methods.

AAAI Conference 2022 Conference Paper

Sublinear Time Approximation of Text Similarity Matrices

• Archan Ray
• Nicholas Monath
• Andrew McCallum
• Cameron Musco

We study algorithms for approximating pairwise similarity matrices that arise in natural language processing. Generally, computing a similarity matrix for n data points requires Ω(n2 ) similarity computations. This quadratic scaling is a significant bottleneck, especially when similarities are computed via expensive functions, e. g. , via transformer models. Approximation methods reduce this quadratic complexity, often by using a small subset of exactly computed similarities to approximate the remainder of the complete pairwise similarity matrix. Significant work focuses on the efficient approximation of positive semidefinite (PSD) similarity matrices, which arise e. g. , in kernel methods. However, much less is understood about indefinite (non-PSD) similarity matrices, which often arise in NLP. Motivated by the observation that many of these matrices are still somewhat close to PSD, we introduce a generalization of the popular Nyström method to the indefinite setting. Our algorithm can be applied to any similarity matrix and runs in sublinear time in the size of the matrix, producing a rank-s approximation with just O(ns) similarity computations. We show that our method, along with a simple variant of CUR decomposition, performs very well in approximating a variety of similarity matrices arising in NLP tasks. We demonstrate high accuracy of the approximated similarity matrices in tasks of document classification, sentence similarity, and cross-document coreference.

NeurIPS Conference 2021 Conference Paper

Coresets for Classification – Simplified and Strengthened

• Tung Mai
• Cameron Musco
• Anup Rao

We give relative error coresets for training linear classifiers with a broad class of loss functions, including the logistic loss and hinge loss. Our construction achieves $(1\pm \epsilon)$ relative error with $\tilde O(d \cdot \mu_y(X)^2/\epsilon^2)$ points, where $\mu_y(X)$ is a natural complexity measure of the data matrix $X \in \mathbb{R}^{n \times d}$ and label vector $y \in \{-1, 1\}^n$, introduced by Munteanu et al. 2018. Our result is based on subsampling data points with probabilities proportional to their $\ell_1$ $Lewis$ $weights$. It significantly improves on existing theoretical bounds and performs well in practice, outperforming uniform subsampling along with other importance sampling methods. Our sampling distribution does not depend on the labels, so can be used for active learning. It also does not depend on the specific loss function, so a single coreset can be used in multiple training scenarios.

ICML Conference 2021 Conference Paper

DeepWalking Backwards: From Embeddings Back to Graphs

• Sudhanshu Chanpuriya
• Cameron Musco
• Konstantinos Sotiropoulos
• Charalampos E. Tsourakakis

Low-dimensional node embeddings play a key role in analyzing graph datasets. However, little work studies exactly what information is encoded by popular embedding methods, and how this information correlates with performance in downstream learning tasks. We tackle this question by studying whether embeddings can be inverted to (approximately) recover the graph used to generate them. Focusing on a variant of the popular DeepWalk method \cite{PerozziAl-RfouSkiena: 2014, QiuDongMa: 2018}, we present algorithms for accurate embedding inversion – i. e. , from the low-dimensional embedding of a graph $G$, we can find a graph $\tilde G$ with a very similar embedding. We perform numerous experiments on real-world networks, observing that significant information about $G$, such as specific edges and bulk properties like triangle density, is often lost in $\tilde G$. However, community structure is often preserved or even enhanced. Our findings are a step towards a more rigorous understanding of exactly what information embeddings encode about the input graph, and why this information is useful for learning tasks.

ICML Conference 2021 Conference Paper

Faster Kernel Matrix Algebra via Density Estimation

• Arturs Backurs
• Piotr Indyk
• Cameron Musco
• Tal Wagner

We study fast algorithms for computing basic properties of an n x n positive semidefinite kernel matrix K corresponding to n points x_1, .. ., x_n in R^d. In particular, we consider the estimating the sum of kernel matrix entries, along with its top eigenvalue and eigenvector. These are some of the most basic problems defined over kernel matrices. We show that the sum of matrix entries can be estimated up to a multiplicative factor of 1+\epsilon in time sublinear in n and linear in d for many popular kernel functions, including the Gaussian, exponential, and rational quadratic kernels. For these kernels, we also show that the top eigenvalue (and a witnessing approximate eigenvector) can be approximated to a multiplicative factor of 1+\epsilon in time sub-quadratic in n and linear in d. Our algorithms represent significant advances in the best known runtimes for these problems. They leverage the positive definiteness of the kernel matrix, along with a recent line of work on efficient kernel density estimation.

NeurIPS Conference 2021 Conference Paper

On the Power of Edge Independent Graph Models

• Sudhanshu Chanpuriya
• Cameron Musco
• Konstantinos Sotiropoulos
• Charalampos Tsourakakis

Why do many modern neural-network-based graph generative models fail to reproduce typical real-world network characteristics, such as high triangle density? In this work we study the limitations of $edge\ independent\ random\ graph\ models$, in which each edge is added to the graph independently with some probability. Such models include both the classic Erdos-Renyi and stochastic block models, as well as modern generative models such as NetGAN, variational graph autoencoders, and CELL. We prove that subject to a $bounded\ overlap$ condition, which ensures that the model does not simply memorize a single graph, edge independent models are inherently limited in their ability to generate graphs with high triangle and other subgraph densities. Notably, such high densities are known to appear in real-world social networks and other graphs. We complement our negative results with a simple generative model that balances overlap and accuracy, performing comparably to more complex models in reconstructing many graph statistics.

ICML Conference 2020 Conference Paper

Efficient Intervention Design for Causal Discovery with Latents

• Raghavendra Addanki
• Shiva Prasad Kasiviswanathan
• Andrew McGregor 0001
• Cameron Musco

We consider recovering a causal graph in presence of latent variables, where we seek to minimize the cost of interventions used in the recovery process. We consider two intervention cost models: (1) a linear cost model where the cost of an intervention on a subset of variables has a linear form, and (2) an identity cost model where the cost of an intervention is the same, regardless of what variables it is on, i. e. , the goal is just to minimize the number of interventions. Under the linear cost model, we give an algorithm to identify the ancestral relations of the underlying causal graph, achieving within a $2$-factor of the optimal intervention cost. This approximation factor can be improved to $1+\eps$ for any $\eps > 0$ under some mild restrictions. Under the identity cost model, we bound the number of interventions needed to recover the entire causal graph, including the latent variables, using a parameterization of the causal graph through a special type of colliders. In particular, we introduce the notion of $p$-colliders, that are colliders between pair of nodes arising from a specific type of conditioning in the causal graph, and provide an upper bound on the number of interventions as a function of the maximum number of $p$-colliders between any two nodes in the causal graph.

SODA Conference 2020 Conference Paper

Fast and Space Efficient Spectral Sparsification in Dynamic Streams

• Michael Kapralov
• Aida Mousavifar
• Cameron Musco
• Christopher Musco
• Navid Nouri
• Aaron Sidford
• Jakab Tardos

In this paper, we resolve the complexity problem of spectral graph sparcification in dynamic streams up to polylogarithmic factors. Using a linear sketch we design a streaming algorithm that uses Õ ( n ) space, and with high probability, recovers a spectral sparsifier from the sketch in Õ ( n ) time. 1 Prior results either achieved near optimal Õ ( n ) space, but Ω( n 2 ) recovery time [Kapralov et al. ‘14], or ran in o ( n 2 ) time, but used polynomially suboptimal space [Ahn et al ‘13]. Our main technical contribution is a novel method for recovering graph edges with high effective resistance from a linear sketch. We show how to do so in nearly linear time by ‘bucketing’ vertices of the input graph into clusters using a coarse approximation to the graph's effective resistance metric. A second main contribution is a new pseudorandom generator (PRG) for linear sketching algorithms. Constructed from a locally computable randomness extractor, our PRG stretches a seed of Õ (n) random bits polynomially in length with just log O (1) n run-time cost per evaluation. This improves on Nisan's commonly used PRG, which in our setting would require Õ ( n ) time per evaluation. Our faster PRG is essential to simultaneously achieving near optimal space and time complexity.

NeurIPS Conference 2020 Conference Paper

Fourier Sparse Leverage Scores and Approximate Kernel Learning

• Tamas Erdelyi
• Cameron Musco
• Christopher Musco

We prove new explicit upper bounds on the leverage scores of Fourier sparse functions under both the Gaussian and Laplace measures. In particular, we study s-sparse functions of the form $f(x) = \sum_{j=1}^s a_j e^{i \lambda_j x}$ for coefficients $a_j \in C$ and frequencies $\lambda_j \in R$. Bounding Fourier sparse leverage scores under various measures is of pure mathematical interest in approximation theory, and our work extends existing results for the uniform measure [Erd17, CP19a]. Practically, our bounds are motivated by two important applications in machine learning: 1. Kernel Approximation. They yield a new random Fourier features algorithm for approximating Gaussian and Cauchy (rational quadratic) kernel matrices. For low-dimensional data, our method uses a near optimal number of features, and its runtime is polynomial in the *statistical dimension* of the approximated kernel matrix. It is the first ``oblivious sketching method'' with this property for any kernel besides the polynomial kernel, resolving an open question of [AKM+17, AKK+20b]. 2. Active Learning. They can be used as non-uniform sampling distributions for robust active learning when data follows a Gaussian or Laplace distribution. Using the framework of [AKM+19], we provide essentially optimal results for bandlimited and multiband interpolation, and Gaussian process regression. These results generalize existing work that only applies to uniformly distributed data.

FOCS Conference 2020 Conference Paper

Near Optimal Linear Algebra in the Online and Sliding Window Models

• Vladimir Braverman
• Petros Drineas
• Cameron Musco
• Christopher Musco
• Jalaj Upadhyay
• David P. Woodruff
• Samson Zhou

We initiate the study of numerical linear algebra in the sliding window model, where only the most recent W updates in a stream form the underlying data set. Although many existing algorithms in the sliding window model use or borrow elements from the smooth histogram framework (Braverman and Ostrovsky, FOCS 2007), we show that many interesting linear-algebraic problems, including spectral and vector induced matrix norms, generalized regression, and low-rank approximation, are not amenable to this approach in the row-arrival model. To overcome this challenge, we first introduce a unified row-sampling based framework that gives randomized algorithms for spectral approximation, low-rank approximation/projection-cost preservation, and ℓ 1 -subspace embeddings in the sliding window model, which often use nearly optimal space and achieve nearly input sparsity runtime. Our algorithms are based on “reverse online” versions of offline sampling distributions such as (ridge) leverage scores, ℓ 1 sensitivities, and Lewis weights to quantify both the importance and the recency of a row; our structural results on these distributions may be of independent interest for future algorithmic design. Although our techniques initially address numerical linear algebra in the sliding window model, our row-sampling framework rather surprisingly implies connections to the well-studied online model; our structural results also give the first sample optimal (up to lower order terms) online algorithm for low-rank approximation/projection-cost preservation. Using this powerful primitive, we give online algorithms for column/row subset selection and principal component analysis that resolves the main open question of Bhaskara et al. (FOCS 2019). We also give the first online algorithm for ℓ 1 -subspace embeddings. We further formalize the connection between the online model and the sliding window model by introducing an additional unified framework for deterministic algorithms using a merge and reduce paradigm and the concept of online coresets, which we define as a weighted subset of rows of the input matrix that can be used to compute a good approximation to some given function on all of its prefixes. Our sampling based algorithms in the row-arrival online model yield online coresets, giving deterministic algorithms for spectral approximation, low-rank approximation/projection-cost preservation, and ℓ 1 - subspace embeddings in the sliding window model that use nearly optimal space.

NeurIPS Conference 2020 Conference Paper

Node Embeddings and Exact Low-Rank Representations of Complex Networks

• Sudhanshu Chanpuriya
• Cameron Musco
• Konstantinos Sotiropoulos
• Charalampos Tsourakakis

Low-dimensional embeddings, from classical spectral embeddings to modern neural-net-inspired methods, are a cornerstone in the modeling and analysis of complex networks. Recent work by Seshadhri et al. (PNAS 2020) suggests that such embeddings cannot capture local structure arising in complex networks. In particular, they show that any network generated from a natural low-dimensional model cannot be both sparse and have high triangle density (high clustering coefficient), two hallmark properties of many real-world networks. In this work we show that the results of Seshadhri et al. are intimately connected to the model they use rather than the low-dimensional structure of complex networks. Specifically, we prove that a minor relaxation of their model can generate sparse graphs with high triangle density. Surprisingly, we show that this same model leads to exact low-dimensional factorizations of many real-world networks. We give a simple algorithm based on logistic principal component analysis (LPCA) that succeeds in finding such exact embeddings. Finally, we perform a large number of experiments that verify the ability of very low-dimensional embeddings to capture local structure in real-world networks.

SODA Conference 2020 Conference Paper

Sample Efficient Toeplitz Covariance Estimation

• Yonina C. Eldar
• Jerry Li 0001
• Cameron Musco
• Christopher Musco

STOC Conference 2019 Conference Paper

A universal sampling method for reconstructing signals with simple Fourier transforms

• Haim Avron
• Michael Kapralov
• Cameron Musco
• Christopher Musco
• Ameya Velingker
• Amir Zandieh

Reconstructing continuous signals based on a small number of discrete samples is a fundamental problem across science and engineering. We are often interested in signals with "simple'' Fourier structure -- e.g., those involving frequencies within a bounded range, a small number of frequencies, or a few blocks of frequencies -- i.e., bandlimited, sparse, and multiband signals, respectively. More broadly, any prior knowledge on a signal's Fourier power spectrum can constrain its complexity. Intuitively, signals with more highly constrained Fourier structure require fewer samples to reconstruct.

NeurIPS Conference 2019 Conference Paper

Toward a Characterization of Loss Functions for Distribution Learning

• Nika Haghtalab
• Cameron Musco
• Bo Waggoner

In this work we study loss functions for learning and evaluating probability distributions over large discrete domains. Unlike classification or regression where a wide variety of loss functions are used, in the distribution learning and density estimation literature, very few losses outside the dominant \emph{log loss} are applied. We aim to understand this fact, taking an axiomatic approach to the design of loss functions for distributions. We start by proposing a set of desirable criteria that any good loss function should satisfy. Intuitively, these criteria require that the loss function faithfully evaluates a candidate distribution, both in expectation and when estimated on a few samples. Interestingly, we observe that \emph{no loss function} possesses all of these criteria. However, one can circumvent this issue by introducing a natural restriction on the set of candidate distributions. Specifically, we require that candidates are \emph{calibrated} with respect to the target distribution, i. e. , they may contain less information than the target but otherwise do not significantly distort the truth. We show that, after restricting to this set of distributions, the log loss and a large variety of other losses satisfy the desired criteria. These results pave the way for future investigations of distribution learning that look beyond the log loss, choosing a loss function based on application or domain need.

NeurIPS Conference 2018 Conference Paper

Inferring Networks From Random Walk-Based Node Similarities

• Jeremy Hoskins
• Cameron Musco
• Christopher Musco
• Babis Tsourakakis

Digital presence in the world of online social media entails significant privacy risks. In this work we consider a privacy threat to a social network in which an attacker has access to a subset of random walk-based node similarities, such as effective resistances (i. e. , commute times) or personalized PageRank scores. Using these similarities, the attacker seeks to infer as much information as possible about the network, including unknown pairwise node similarities and edges. For the effective resistance metric, we show that with just a small subset of measurements, one can learn a large fraction of edges in a social network. We also show that it is possible to learn a graph which accurately matches the underlying network on all other effective resistances. This second observation is interesting from a data mining perspective, since it can be expensive to compute all effective resistances or other random walk-based similarities. As an alternative, our graphs learned from just a subset of effective resistances can be used as surrogates in a range of applications that use effective resistances to probe graph structure, including for graph clustering, node centrality evaluation, and anomaly detection. We obtain our results by formalizing the graph learning objective mathematically, using two optimization problems. One formulation is convex and can be solved provably in polynomial time. The other is not, but we solve it efficiently with projected gradient and coordinate descent. We demonstrate the effectiveness of these methods on a number of social networks obtained from Facebook. We also discuss how our methods can be generalized to other random walk-based similarities, such as personalized PageRank scores. Our code is available at https: //github. com/cnmusco/graph-similarity-learning.

SODA Conference 2018 Conference Paper

Stability of the Lanczos Method for Matrix Function Approximation

• Cameron Musco
• Christopher Musco
• Aaron Sidford

Theoretically elegant and ubiquitous in practice, the Lanczos method can approximate f ( A ) x for any symmetric matrix A ∊ ℝ n × n, vector x ∊ ℝ n, and function f. In exact arithmetic, the method's error after k iterations is bounded by the error of the best degree- k polynomial uniformly approximating the scalar function f ( x ) on the range [λ min ( A ), λ max ( A )]. However, despite decades of work, it has been unclear if this powerful guarantee holds in finite precision. We resolve this problem, proving that when, Lanczos essentially matches the exact arithmetic guarantee if computations use roughly log( nC ║ A ║) bits of precision. Our proof extends work of Druskin and Knizhnerman [11], leveraging the stability of the classic Chebyshev recurrence to bound the stability of any polynomial approximating f ( x ). We also study the special case of f ( A ) = A –1 for positive definite A, where stronger guarantees hold for Lanczos. In exact arithmetic the algorithm performs as well as the best polynomial approximating 1/ x at each of A 's eigenvalues, rather than on the full range [λ min ( A ), λ max ( A )]. In seminal work, Greenbaum gives a natural approach to extending this bound to finite precision: she proves that finite precision Lanczos and the related conjugate gradient method match any polynomial approximating 1/ x in a tiny range around each eigenvalue [17]. For A –1, Greenbaum's bound appears stronger than our result. However, we exhibit matrices with condition number κ where exact arithmetic Lanczos converges in polylog( k ) iterations, but Greenbaum's bound predicts at best Ω( κ 1/5 ) iterations in finite precision. It thus cannot offer more than a polynomial improvement over the O ( κ 1/2 ) bound achievable via our result for general f ( A ). Our analysis bounds the power of stable approximating polynomials and raises the question of if they fully characterize the behavior of finite precision Lanczos in solving linear systems. If they do, convergence in less than poly( K ) iterations cannot be expected, even for matrices with clustered, skewed, or otherwise favorable eigenvalue distributions.

SODA Conference 2017 Conference Paper

Input Sparsity Time Low-rank Approximation via Ridge Leverage Score Sampling

• Michael B. Cohen
• Cameron Musco
• Christopher Musco

We present a new algorithm for finding a near optimal low-rank approximation of a matrix A in O ( n nz (A)) time. Our method is based on a recursive sampling scheme for computing a representative subset of A ‘s columns, which is then used to find a low-rank approximation. This approach differs substantially from prior O ( n nz (A)) time algorithms, which are all based on fast Johnson-Lindenstrauss random projections. Our algorithm matches the guarantees of the random projection methods while offering a number of advantages. In addition to better performance on sparse and structured data, sampling algorithms can be applied in settings where random projections cannot. For example, we give new streaming algorithms for the column subset selection and projection-cost preserving sample problems. Our method has also been used in the fastest algorithms for provably accurate Nyström approximation of kernel matrices [56].

NeurIPS Conference 2017 Conference Paper

Is Input Sparsity Time Possible for Kernel Low-Rank Approximation?

• Cameron Musco
• David Woodruff

Low-rank approximation is a common tool used to accelerate kernel methods: the $n \times n$ kernel matrix $K$ is approximated via a rank-$k$ matrix $\tilde K$ which can be stored in much less space and processed more quickly. In this work we study the limits of computationally efficient low-rank kernel approximation. We show that for a broad class of kernels, including the popular Gaussian and polynomial kernels, computing a relative error $k$-rank approximation to $K$ is at least as difficult as multiplying the input data matrix $A \in R^{n \times d}$ by an arbitrary matrix $C \in R^{d \times k}$. Barring a breakthrough in fast matrix multiplication, when $k$ is not too large, this requires $\Omega(nnz(A)k)$ time where $nnz(A)$ is the number of non-zeros in $A$. This lower bound matches, in many parameter regimes, recent work on subquadratic time algorithms for low-rank approximation of general kernels [MM16, MW17], demonstrating that these algorithms are unlikely to be significantly improved, in particular to $O(nnz(A))$ input sparsity runtimes. At the same time there is hope: we show for the first time that $O(nnz(A))$ time approximation is possible for general radial basis function kernels (e. g. , the Gaussian kernel) for the closely related problem of low-rank approximation of the kernelized dataset.

ICML Conference 2017 Conference Paper

Random Fourier Features for Kernel Ridge Regression: Approximation Bounds and Statistical Guarantees

• Haim Avron
• Michael Kapralov
• Cameron Musco
• Christopher Musco
• Ameya Velingker
• Amir Zandieh

Random Fourier features is one of the most popular techniques for scaling up kernel methods, such as kernel ridge regression. However, despite impressive empirical results, the statistical properties of random Fourier features are still not well understood. In this paper we take steps toward filling this gap. Specifically, we approach random Fourier features from a spectral matrix approximation point of view, give tight bounds on the number of Fourier features required to achieve a spectral approximation, and show how spectral matrix approximation bounds imply statistical guarantees for kernel ridge regression.

NeurIPS Conference 2017 Conference Paper

Recursive Sampling for the Nystrom Method

• Cameron Musco
• Christopher Musco

We give the first algorithm for kernel Nystrom approximation that runs in linear time in the number of training points and is provably accurate for all kernel matrices, without dependence on regularity or incoherence conditions. The algorithm projects the kernel onto a set of s landmark points sampled by their ridge leverage scores, requiring just O(ns) kernel evaluations and O(ns^2) additional runtime. While leverage score sampling has long been known to give strong theoretical guarantees for Nystrom approximation, by employing a fast recursive sampling scheme, our algorithm is the first to make the approach scalable. Empirically we show that it finds more accurate kernel approximations in less time than popular techniques such as classic Nystrom approximation and the random Fourier features method.

FOCS Conference 2017 Conference Paper

Sublinear Time Low-Rank Approximation of Positive Semidefinite Matrices

• Cameron Musco
• David P. Woodruff

We show how to compute a relative-error low-rank approximation to any positive semidefinite (PSD) matrix in sublinear time, i. e. , for any n x n PSD matrix A, in Õ(n · poly(k/ε)) time we output a rank-k matrix B, in factored form, for which ||A - B|| 2 F ≤ (1 + ε)||A - A k || F 2, where Ak is the best rank-k approximation to A. When k and 1/ε are not too large compared to the sparsity of A, our algorithm does not need to read all entries of the matrix. Hence, we significantly improve upon previous nnz(A) time algorithms based on oblivious subspace embeddings, and bypass an nnz(A) time lower bound for general matrices (where nnz(A) denotes the number of non-zero entries in the matrix). We prove time lower bounds for low-rank approximation of PSD matrices, showing that our algorithm is close to optimal. Finally, we extend our techniques to give sublinear time algorithms for lowrank approximation of A in the (often stronger) spectral norm metric ||A - B|| 2 2 and for ridge regression on PSD matrices.

ICML Conference 2016 Conference Paper

Faster Eigenvector Computation via Shift-and-Invert Preconditioning

• Dan Garber
• Elad Hazan
• Chi Jin 0001
• Sham M. Kakade
• Cameron Musco
• Praneeth Netrapalli
• Aaron Sidford

We give faster algorithms and improved sample complexities for the fundamental problem of estimating the top eigenvector. Given an explicit matrix $A \in \mathbb{R}^{n \times d}$, we show how to compute an $\epsilon$-approximate top eigenvector of $A^TA$ in time $\tilde O\left( \left[\text{nnz}(A) + \frac{d \text{sr}(A)}{\text{gap}^2} \right] \cdot \log 1/\epsilon\right)$. Here $\text{nnz}(A)$ is the number of nonzeros in $A$, $\text{sr}(A)$ is the stable rank, and gap is the relative eigengap. We also consider an online setting in which, given a stream of i. i. d. samples from a distribution D with covariance matrix $\Sigma$ and a vector $x_0$ which is an $O(\text{gap})$ approximate top eigenvector for $\Sigma$, we show how to refine $x_0$ to an $\epsilon$ approximation using $O \left( \frac{\text{var}(\mathcal{D})}{\text{gap}-\epsilon}\right)$ samples from $\mathcal{D}$. Here $\text{var}(\mathcal{D})$ is a natural notion of variance. Combining our algorithm with previous work to initialize $x_0$, we obtain improved sample complexities and runtimes under a variety of assumptions on D. We achieve our results via a robust analysis of the classic shift-and-invert preconditioning method. This technique lets us reduce eigenvector computation to approximately solving a series of linear systems with fast stochastic gradient methods.

ICML Conference 2016 Conference Paper

Principal Component Projection Without Principal Component Analysis

• Roy Frostig
• Cameron Musco
• Christopher Musco
• Aaron Sidford

We show how to efficiently project a vector onto the top principal components of a matrix, *without explicitly computing these components*. Specifically, we introduce an iterative algorithm that provably computes the projection using few calls to any black-box routine for ridge regression. By avoiding explicit principal component analysis (PCA), our algorithm is the first with no runtime dependence on the number of top principal components. We show that it can be used to give a fast iterative method for the popular principal component regression problem, giving the first major runtime improvement over the naive method of combining PCA with regression. To achieve our results, we first observe that ridge regression can be used to obtain a "smooth projection" onto the top principal components. We then sharpen this approximation to true projection using a low-degree polynomial approximation to the matrix step function. Step function approximation is a topic of long-term interest in scientific computing. We extend prior theory by constructing polynomials with simple iterative structure and rigorously analyzing their behavior under limited precision.

STOC Conference 2015 Conference Paper

Dimensionality Reduction for k-Means Clustering and Low Rank Approximation

• Michael B. Cohen
• Sam Elder
• Cameron Musco
• Christopher Musco
• Madalina Persu

We show how to approximate a data matrix A with a much smaller sketch ~A that can be used to solve a general class of constrained k-rank approximation problems to within (1+ε) error. Importantly, this class includes k-means clustering and unconstrained low rank approximation (i.e. principal component analysis). By reducing data points to just O(k) dimensions, we generically accelerate any exact, approximate, or heuristic algorithm for these ubiquitous problems. For k-means dimensionality reduction, we provide (1+ε) relative error results for many common sketching techniques, including random row projection, column selection, and approximate SVD. For approximate principal component analysis, we give a simple alternative to known algorithms that has applications in the streaming setting. Additionally, we extend recent work on column-based matrix reconstruction, giving column subsets that not only 'cover' a good subspace for A}, but can be used directly to compute this subspace. Finally, for k-means clustering, we show how to achieve a (9+ε) approximation by Johnson-Lindenstrauss projecting data to just O(log k/ε 2 ) dimensions. This is the first result that leverages the specific structure of k-means to achieve dimension independent of input size and sublinear in k.

NeurIPS Conference 2015 Conference Paper

Randomized Block Krylov Methods for Stronger and Faster Approximate Singular Value Decomposition

• Cameron Musco
• Christopher Musco

Since being analyzed by Rokhlin, Szlam, and Tygert and popularized by Halko, Martinsson, and Tropp, randomized Simultaneous Power Iteration has become the method of choice for approximate singular value decomposition. It is more accurate than simpler sketching algorithms, yet still converges quickly for any matrix, independently of singular value gaps. After ~O(1/epsilon) iterations, it gives a low-rank approximation within (1+epsilon) of optimal for spectral norm error. We give the first provable runtime improvement on Simultaneous Iteration: a randomized block Krylov method, closely related to the classic Block Lanczos algorithm, gives the same guarantees in just ~O(1/sqrt(epsilon)) iterations and performs substantially better experimentally. Our analysis is the first of a Krylov subspace method that does not depend on singular value gaps, which are unreliable in practice. Furthermore, while it is a simple accuracy benchmark, even (1+epsilon) error for spectral norm low-rank approximation does not imply that an algorithm returns high quality principal components, a major issue for data applications. We address this problem for the first time by showing that both Block Krylov Iteration and Simultaneous Iteration give nearly optimal PCA for any matrix. This result further justifies their strength over non-iterative sketching methods.

FOCS Conference 2014 Conference Paper

Single Pass Spectral Sparsification in Dynamic Streams

• Michael Kapralov
• Yin Tat Lee
• Cameron Musco
• Christopher Musco
• Aaron Sidford

We present the first single pass algorithm for computing spectral sparsifiers of graphs in the dynamic semi-streaming model. Given a single pass over a stream containing insertions and deletions of edges to a graph, G, our algorithm maintains a randomized linear sketch of the incidence matrix into dimension O(1/∈ 2n polylog(n)). Using this sketch, the algorithm can output a (1±∈) spectral sparsifier for G with high probability. While O(1/∈ 2n polylog(n)) space algorithms are known for computing cut sparsifiers in dynamic streams [1], [2] and spectral sparsifiers in insertion-only streams [3], prior to our work, the best known single pass algorithm for maintaining spectral sparsifiers in dynamic streams required sketches of dimension Ω(1/∈ 2n5/3 ). To achieve our result, we show that, using a coarse sparsifier of G and a linear sketch of G's incidence matrix, it is possible to sample edges by effective resistance, obtaining a spectral sparsifier of arbitrary precision. Sampling from the sketch requires a novel application of ℓ 2 /ℓ 2 sparse recovery, a natural extension of the ℓ 0 methods used for cut sparsifiers in [1]. Recent work of [2] on row sampling for matrix approximation gives a recursive approach for obtaining the required coarse sparsifiers. Under certain restrictions, our approach also extends to the problem of maintaining a spectral approximation for a general matrix A T A given a stream of updates to rows in A.