On Solving Linear Systems in Sublinear Time Alexandr Andoni, - PowerPoint PPT Presentation

On Solving Linear Systems in Sublinear Time Alexandr Andoni, Columbia University Robert Krauthgamer, Weizmann Institute Yosef Pogrow, Weizmann Institute  Google WOLA 2019

Solving Linear Systems Input: 𝐵 ∈ ℝ 𝑜×𝑜 and 𝑐 ∈ ℝ 𝑜  Output: vector 𝑦 that solves 𝐵𝑦 = 𝑐  Many algorithms, different variants:  Matrix 𝐵 is sparse, Laplacian, PSD etc.  Bounded precision (solution 𝑦 is approximate) vs. exact arithmetic  Significant progress: Linear system in Laplacian matrix 𝑀 𝐻 can be  1 solved approximately in near-linear time ෨ 𝑃(nnz 𝑀 𝐻 ⋅ log 𝜗 ) [Spielman- Teng ’ 04, … , Cohen-Kyng-Miller-Pachocky-Peng-Rao-Xu ’ 14] Our focus: Sublinear running time On Solving Linear Systems in Sublinear Time 2

Sublinear-Time Solver Input: 𝐵 ∈ ℝ 𝑜×𝑜 , 𝑐 ∈ ℝ 𝑜 (also 𝜗 > 0 ) and 𝑗 ∈ [𝑜]  𝑦 𝑗 from (any) solution 𝑦 ∗ to 𝐵𝑦 = 𝑐 Output: approximate coordinate ො  𝑦 − 𝑦 ∗ ∞ ≤ 𝜗 𝑦 ∗ ∞ Accuracy bound ො  Formal requirement: There is a solution 𝑦 ∗ to the system, such that  ∗ ≤ 𝜗 𝑦 ∗ ∞ ≥ 3 ∀𝑗 ∈ 𝑜 , Pr ො 𝑦 𝑗 − 𝑦 𝑗 4 Follows framework of Local Computation Algorithms (LCA),  previously used for graph problems [Rubinfeld-Tamir-Vardi-Xie ’ 10] On Solving Linear Systems in Sublinear Time 3

Motivation Fast quantum algorithms for solving linear systems and for machine  learning problems [Harrow-Hassidim-Lloyd ’ 09, … ] Can we match their performance classically?  Recent success story: quantum  classical algorithm [Tang ’ 18]  New direction in sublinear-time algorithms  “ Local ” computation in numerical problems  Compare computational models (representation, preprocessing),  accuracy guarantees, input families (e.g., Laplacian vs. PSD) Known quantum algorithms have modeling requirements (e.g., quantum  encoding of 𝑐 ) On Solving Linear Systems in Sublinear Time 4

Algorithm for Laplacians Informally: Can solve Laplacian systems of bounded-degree  expander in polylog(n) time Key limitations: sparsity and condition number  Notation:  𝑀 𝐻 = 𝐸 − 𝐵 is the Laplacian matrix of graph 𝐻  + is its Moore-Penrose pseudo-inverse 𝑀 𝐻  Theorem 1: Suppose the input is a 𝑒 -regular 𝑜 -vertex graph 𝐻 ,  together with its condition number 𝜆 > 0 , 𝑐 ∈ ℝ 𝑜 , 𝑣 ∈ 𝑜 and 𝜗 > 0 . 𝑦 𝑣 ∈ ℝ such that for 𝑦 ∗ = 𝑀 𝐻 + 𝑐 , Our algorithm computes ො ∗ ≤ 𝜗 𝑦 ∗ ∞ ≥ 3 ∀𝑣 ∈ 𝑜 , Pr ො 𝑦 𝑣 − 𝑦 𝑣 4 , and runs in time ෨ 𝑃(𝑒𝜗 −2 𝑡 3 ) for 𝑡 = ෨ 𝑃(𝜆 log 𝑜) . More inputs? Faster? On Solving Linear Systems in Sublinear Time 5

Some Extensions Can replace 𝑜 with 𝑐 0  Example: Effective resistance can be approximate (in expanders) in  constant running time! 𝑆eff(𝑣, 𝑤) = 𝑓 𝑣 − 𝑓 𝑤 𝑈 𝑀 𝐻 + (𝑓 𝑣 − 𝑓 𝑤 ) Improved running time if  Graph 𝐻 is preprocessed  One can sample a neighbor in 𝐻 , or  Extends to Symmetric Diagonally Dominant (SDD) matrix 𝑇  𝜆 is condition number of 𝐸 −1/2 𝑇𝐸 −1/2  On Solving Linear Systems in Sublinear Time 6

Lower Bound for PSD Systems Informally: Solving “ similar ” PSD systems requires polynomial time  Similar = bounded condition number and sparsity  Even if the matrix can be preprocessed  Theorem 2: For certain invertible PSD matrices 𝑇 , with bounded  sparsity 𝑒 and condition number 𝜆 , every randomized algorithm must query 𝑜 Ω(1/𝑒 2 ) coordinates of the input 𝑐 . Here, the output is ො 𝑦 𝑣 ∈ ℝ for a fixed 𝑣 ∈ 𝑜 , required to satisfy  ∗ ≤ 1 3 5 𝑦 ∗ ∞ ≥ ∀𝑣 ∈ 𝑜 , Pr 𝑦 𝑣 − 𝑦 𝑣 ො 4 , for 𝑦 ∗ = 𝑇 −1 𝑐 . In particular, 𝑇 may be preprocessed  On Solving Linear Systems in Sublinear Time 7

Dependence on Condition Number Informally: Quadratic dependence on 𝜆 is necessary  Our algorithmic bound ෩ O(𝜆 3 ) is near-optimal, esp. when matrix 𝑇 can be  preprocessed Theorem 3: For certain graphs 𝐻 of maximum degree 4 and any  condition number 𝜆 > 0 , every randomized algorithm (for 𝑀 𝐻 ) with 1 log 𝑜 must probe ෩ Ω(𝜆 2 ) coordinates of the input 𝑐 . accuracy 𝜗 = Again, the output is ො 𝑦 𝑣 ∈ ℝ for a fixed 𝑣 ∈ 𝑜 , required to satisfy  ∗ ≤ 1 3 log 𝑜 𝑦 ∗ ∞ ≥ ∀𝑣 ∈ 𝑜 , Pr 𝑦 𝑣 − 𝑦 𝑣 ො 4 , for 𝑦 ∗ = 𝑀 𝐻 + 𝑐 . In particular, 𝐻 may be preprocessed  On Solving Linear Systems in Sublinear Time 8

Algorithmic Techniques Famous Monte-Carlo method of von Neumann and Ulam:  Write matrix inverse by power series 𝐽 − 𝑌 −1 = σ 𝑢≥0 𝑌 𝑢 ∀ 𝑌 < 1 , then estimate it by random walks (in 𝑌 ) with unbiased expectation Inverting a Laplacian 𝑀 𝐻 = 𝑒𝐽 − 𝐵 corresponds to summing walks in 𝐻  𝑈 σ 𝑢≥0 𝐵 𝑢 𝑐 as sum over all walks, estimate it by sampling For us: view 𝑓 𝑣  (random walks) Need to control: number of walks and their length  Large powers 𝑢 > 𝑢 ∗ contribute relatively little (by condition number)  Estimate truncated series ( 𝑢 ≤ 𝑢 ∗ ) by short random walks (by Chebyshev ’ s  inequality) On Solving Linear Systems in Sublinear Time 9

Related Work – All Algorithmic Similar techniques were used before in related contexts but under  different assumptions, models and analyses: + [Doron-Le Gall- Probabilistic log-space algorithms for approximating 𝑀 𝐻  Ta-Shma ’ 17] Asks for entire matrix, uses many long random walks (independent of 𝜆 )  Local solver for Laplacian systems with boundary conditions [Chung-  Simpson ’ 15] Solver relies on a different power series and random walks  Local solver for PSD systems [Shyamkumar-Banerjee-Lofgren ’ 16]  Polynomial time nnz 𝑇 2/3 under assumptions like bounded matrix norm and  random 𝑣 ∈ 𝑜 Local solver for Pagerank [Bressan-Peserico-Pretto ’ 18, Borgs-Brautbar-  Chayes-Teng ’ 14] Polynomial time O(𝑜 2/3 ) and O( nd 1/2 ) for certain matrices (non-symmetric  but by definition are diagonally-dominant) On Solving Linear Systems in Sublinear Time 10

Lower Bound Techniques PSD lower bound: Take Laplacian of 2𝑒 -regular expander but with:  𝑣 high girth,  edges signed ±1 at random, and  𝑃( 𝑒) on the diagonal (PSD but not Laplacian)  𝑠 The graph looks like a tree locally  Up to radius Θ log 𝑜 around 𝑣  𝑐 𝑥 = 0 Set 𝑐 𝑥 = ±1 for 𝑥 at distance 𝑠 , and 0 otherwise  𝑐 𝑥 = ±1 Signs have small bias 𝜀 ≈ 𝑒 −𝑠/2  Recovering it requires reading Ω(𝜀 −2 ) entries  Using inversion formula, 𝑦 𝑣 ≈ average of 𝑐 𝑥 ‘ s 𝑐 𝑥 = 0  Condition number lower bound: Take two 3-regular expanders  connected by a matching of size 𝑜/𝜆 Let 𝑐 𝑥 = ±1 with slight bias inside each expander  On Solving Linear Systems in Sublinear Time 11

Further Questions Accuracy guarantees  Different norms?  Condition number of 𝑇 instead of 𝐸 −1/2 𝑇𝐸 −1/2 ?  Other representations (input/output models)?  Access the input 𝑐 via random sampling?  Sample from the output 𝑦 ?  Other numerical problems?  Thank You! On Solving Linear Systems in Sublinear Time 12