Fast Tensor Completion via Approximate Richardson Iteration
Authors: Mehrdad Ghadiri
Matthew Fahrbach
Yunbum Kook
Ali Jadbabaie
Venue: ICML 2025 poster (2025)
Keywords: Tensor completion
tensor decomposition
leverage score sampling
sketching
iterative methods
linear regression
Primary area: optimization
TL;DR: A new tensor completion method that combines leverage score sampling and iterative methods with provable guarantees.

Abstract: We study tensor completion (TC) through the lens of low-rank tensor decomposition (TD). Many TD algorithms use fast alternating minimization methods to solve _highly structured_ linear regression problems at each step (e.g., for CP, Tucker, and tensor-train decompositions). However, such algebraic structure is often lost in TC regression problems, making direct extensions unclear.  This work proposes a novel _lifting_ method for approximately solving TC regression problems using structured TD regression algorithms as blackbox subroutines, enabling sublinear-time methods. We analyze the convergence rate of our approximate Richardson iteration-based algorithm, and our empirical study shows that it can be 100x faster than direct methods for CP completion on real-world tensors.

Fast Tensor Completion via Approximate Richardson Iteration
Mehrdad Ghadiri'
Matthew Fahrbach 2
Yunbum Kook 3 Ali Jadbabaie
Abstract
regularization term. For brevity, we write X without explicit dependence on O. Rank constraints can be incorporated into (1) by including appropriate penalty terms in R(0), c.g., by using & or l2 regularization terms to reduce the effective dimension (Fahrbach et al., 2022).
We study tensor completion (TC) through the lens of low-rank tensor decomposition (TD). Many TD algorithms use fast alternating minimization meth- ods to solve highly structured linear regression problems at each step (e.g., for CP, Tucker, and tensor-train decompositions). However, such alge- braic structure is often lost in TC regression prob- lems, making direct extensions unclear. This work proposes novel lifting method for approximately solving TC regression problems using structured TD regression algorithms as blackbox subroutines, enabling sublinear-time methods. We analyze the convergence rate of our approximate Richardson iteration- -based algorithm, and our empirical study shows that it can be 100x faster than direct meth- ods for CP completion real-world tensors.
This work focuses on the sum of squared errors
()()) = Σ
(fi...i Fnn.,x) (2)
(iI...
To introduce our approach, we use running example where x 9 X1 A (1) X2 A (2) X A (N) is low-rank Tucker decomposition with core tensor G E RRRxx x--xxxxxx Rv and factor matrices A (m) R f. XR,, We do not focus on regulariza- tion in this paper, but our approach extends to regularized problems, provided there exists an algorithm that accelerates the corresponding structured and regularized objectives. For example, Fahrbach et al. (2022) present subquadratic-time algorithms for Kronecker ridge regression (i.e., & regular- ization).
1. Introduction
Tensor completion (TC) is the higher-order generalization of matrix completion. It has many applications across data mining, machine learning, signal processing, and statistics (see Song et al. (2019) for a detailed survey of applications). In the TC problem, we are given a partially observed tensor as input (i.e., only a subset of its entries known), and the goal is to impute the missing values. Under low-rank and other statistical assumptions, we can recover the missing values by minimizing a loss function over only the observed entries.
In general, TC is a nonconvex and NP-hard problem (Hillar & Lim, 2013). The special case where = [41] x x [Iv] is the widely studied tensor decomposition (TD) problem, for which alternating least squares (ALS) methods often used to compute 9, (1) A (N).
Each ALS step fixes all but one of ( ((N), and optimizes the unfixed component. the case of TD, each ALS step is a highly structured linear regression problem. For example, in the core tensor update for 9. ALS solves
Consider a tensor X E IR.4 x.xx with a subset of observed indices C [T,] x ... x [Ix]. The general TC problem is
(@A
N A (n) voc(S') vec(X)) -o(('')voc(33)))
min S'ERRI XOOK
(3)
min (0()) R(0),
(1)
where is the Kronecker product and vec(S") R is the flattened version tensor S', with R = II R..
where O are the learnable parameters, (0) E R4xx-xxx is the reconstructed tensor, Ln is the loss function defining the error between x()) and X on entries in Ω and R()) is the
Diao al. (2019) and Fahrbach al. (2022) recently ex- ploited the Kronecker product structure in (3) to give algo- rithms with running times that are sublinear in the size of the full tensor 7 = IT... 1 T, These approaches use leverage score sampling to approximately solve (3).
"qqual contribution 'MIT "Google Research GGeorgia Tech. Correspondence to: Mehrdad Ghadiri <mchrdadg@@nn @ii..edu>
Proceedings of the 42nd International Conference Machine Learning. Vancouver, Canada. PMLR 267, 2025. Copyright 2025 by the author(s).
In TC, however, only a subset of observations appear in the loss function. Letting A 8 A (a) and b = veo(I))

---

Fast Tensor Complction via Approximate Richardson Iteration
the core tensor update becomes
In Section 3, we propose using the mini-ALS algorithm for each step of ALS for TC. We show that it simulates the preconditioned Richardson iteration (Lemma 3.5), Our main theoretical contribution is proving that the second step of a mini-ALSS iteration can be performed approximately. We quantify how small the approxima- tion error must be for approximate-minn mini---L to con- verge at the same rate as the Richardson iteration (The- orem 3.7). This lets us extend recent line of work on leverage score sampling-based based TD ALS algorithms to the TC setting.
vec ( < arg min - |Aag--blll.
(4)
XERR
where Ag is a submatrix of A whose rows correspond to the indices in S2 Observe that the Kronecker product structure in (3) no longer exists for A in (4), hence fast TD methods do not immediately extend to Q-masked TC versions.
A natural idea overcome the lack of structure in the Ω masked updates is to lift (4) to a nigher-diiennionaa problem by introducing variables by. where S2 is the complement of Ω2 i.e., letting the unobserved entries in X be free variables. This is a convex problem with much of the same structure as the design matrix in the full TD problem. Further, it gives the same solution as the TC update in (4):
In Section 4, we use state-of-the-arr TD ALS algo- rithms for CP, Tucker, and tensor-train decompositions as blackbox subroutines in Algorithm I to obtain novel sampling- -baeed TC algorithms. Hence, our lifting ap- proach for TC also benefits from future TD algorithmic improvements.
(x", ba²) = arg min| Ax - bll2-
(5)
b,
In Section 5, we show that leverage score sampling is an effective method for solving large structured regres- sion problems via the coupled matrix problem. Then we compare the empirical performance of our lifted al- gorithm to direct methods for low-rank CP completion on synthetic and real-world tensors. We observe that mini-ALS can be orders of magnitude faster than direct ALS methods, while achieving comparable reconstruc- tion errors. Finally, we propose an accelerated version with adaptive step sizes that extrapolates the trajectory of the iterates (k) and converges in fewer iterations.
To our knowledge, such ideas date back to Healy & Westma- cot (1956) in the experimental design and causal inference literature. However, it is not clear a priori that lifting is help- ful for computational efficiency-it it restores the structure of the design matrix, but the new problem (5) is larger and higher dimensional.
This work proposes solving (5) with a two-step procedure called mini-ALS. Given a vector ((x) corresponding to an iterate of a block of variables in the low-rank decomposition, it repeats the following:
1.2. Related Work
1. Set (k)
by <- A App (8) and b <- bu
(dc)
// lift
CP completion. Tomasi & Bro (2005) proposed an ALS algorithm for CP completion that repeats the following two- step process: (1) fill in the missing values using the current CP decomposition, and (2) update one factor matrix. Their algorithm is equivalent to running one iteration of mini-ALS in each step of ALS. As Tomasi & Bro (2005) discuss, this can lead slower convergence and an increased likelihood of converging to suboptimal local minima because of errors introduced by the imputed missing values.
2. Set x(++1) + arg min Axx - |00 // solve
Mini-ALS iterations can still be expensive since A R is tall-and-skinny matrix (i.e., I > R). However, step two is a structured ALS for TD, so we propose solving it approximately with row sampling, allowing us to tap into a rich line of work on leverage score sampling for CP de- composition (Cheng et al., 2016; Larsen & Kolda, 2022; Bharadwaj et al., 2023), Tucker decomposition (Diao et al., 2019; Fahrbach et al., 2022), and tensor-train decomposi- tion (Bharadwaj et al., 2024).
In contrast, approximate-mint! in each step of ALS. By doing this, we establish a connec- tion to the Richardson iteration and build on its convergence guarantees. Further, Tomasi & Bro (2005) explicitly fill in all missing values using the current decomposition, whereas we only impute missing values required by row sampling, allowing us to achieve sublinear-time updates in the size of the tensor. In general, iteratively fitting to imputed missing values falls under the umbrella of expectation expectation-maxiiizzttooo maximization (FM) algorithms (Little & Rubin, 2019, Chapter 8).
runs until convergence
Returning to the running Tucker completion example, up- dating the core tensor in step two of mini-ALS via (3) only requires sampling Ô(R)) rows of A, which is a substantially smaller problem. Further, we can compute (k) lazily, i.c., only the entries corresponding to sampled rows. We call our lifted iterative method approximate-minn-ALS...
1.1. Our Contributions
Statistical assumptions. Similar to minimizing the nu- clear norm for matrix completion (Fazel, 2002; Candes &
We summarize our main contributions as follows:
2

---

Fast Tensor Complction via Approximate Richardson Iteration
Recht, 2012), line of research in noisy TC proposes mini- mizing a convex relaxation of rank and identifies statistical assumptions under which the problem is recoverable (Barak & Moitra, 2016). Two standard assumptions are incoher- ence and the missing-at teanng-at-rannden randoam assumption. In Section 3, we discuss how these assumptions provide a bound on the number of steps required for mini-ALS. Alternating mini- mization approaches have also been applied in this line of noisy TC research (Jain & Oh, 2014; Liu & Moitra, 2020). It is also consistently observed that methods based on TD and tensor unfoldings are more practical and computation- ally efficient (Acar et al., 2011; Montanari & Sun, 2018; Filipović& & Jukié, 2015; Shah & Yu, 2019; 2023).
The inner product of two tensors X, y E IR x is
I
IN
(),
=
It Eiybandk kando
The Frobenius norm of a tensor X is ||0||== V ())
2.1. Tensor Decompositions
The tensor decompositions of X 1 below can be seen as higher-order analogs of low-rank matrix factor- ization. We direct the reader to Kolda & Bader (2009) for comprehensive survey on this topic.
2. Preliminaries
CP decomposition. A rank-R CP decomposition repre- sents X with À E R and N factors A A(r) E R for n E [N), where each column of A (n) has unit norm. The reconstructed tensor X is defined elementwise as:
Notation. The order of tensor is the number of its di- mensions N. We denote scalars by normal lowercase letters E ER vectors by boldface lowercase letters x E RR", matri- ces by boldface uppercase letters X E mm and higher- order tensors by boldface script letters E xx We use normal uppercase letters for the size of an index set, e.g., [N] {1, ,,...,,,, NN. We define 1+ :- I/ for k [N] and similarly R Re- We denote the i-th entry of x by rr, the ji, entry of X by Tij. and the (i, k) entry of a third-order tensor X by Tsjk-
R
(N) Giyr
i AN =
Ap Air
Tucker decomposition. A rank-r Tucker decomposition represents X with a core tensor 9 E IR x--xXx and N factor matrices A (n) E R x R., for n E [N], where r = , Ry) is the multilinear rank (Ghadiri et al., 2023a). The reconstructed tensor X = 9 x1 A X2 ... XN A (N) is defined elementwise as:
Linear algebra. A symmetric matrix A E Rax= is pos- itive semi-definite (PSD) if v Av 2 for any v E R For two symmetric matrices A, B, we use & B in- dicate that B- A is PSD. For PSD matrix M E Ruxn and vector v E RR, we define| | M == (v Mv) 1/2 For A E R b E R and f C [m]. we use Aq E R [0|xn and bn E R [ to denote the submatrix and subvector with rows indexed by Ω. We let @ denote the Kronecker product and @ denote the Khatri- Rao product.
R
Ry gr a
(1)
(N)
Biy..dy =
-.
iiry
ry=1
TT decomposition. A rank-r tensor train (TT) decompo- sition (Oseledets, 2011) represents X with N third-order TT-cores A (N) E IR PRa-1xlax Rs, for n E [N], using the convention Ro Ry 1. The reconstructed tensor X is defined elementwise as:
Tensor algebra. The fibers of tensor are the vectors ob- tained by fixing all but one index, e.g., if X E X the column, row and tube fibers are X.jk. Xit and X;j:, respec- tively. The mode-n unfolding of a tensor X E R is the matrix X(n) E R (II 1/a 1>) that arranges the mode-n fibers of X as rows of X (m) sorted lexicograph- ically by index. The vectorization vec (X) E h o stacks the elements of X lexicographically by index.
Big..oi
(N-1) カNー
(N)
ANN
1x R R2
R
Ry-1xt
Remark 2.1. All three of these TDs are instances of the more general tensor network framework (see Appendix A.4).
For n. E [N], we denote the mode-n product of tensor X E R/ XIx and matrix A E J.I by y x xn A, where y E R4x--x1, 1xJxIn (1 xx-xIx This operation multiplies each mode-n fiber of by A, and is expressed elementwise as
2.2. ALS Formulations
Alternating least squares (ALS) methods are the gold stan- dard for low-rank tensor decomposition, e.g., they are the first techniques mentioned in the MATLAB Tensor Toolbox (Bader & Kolda, 2023). ALS cyclically minimizes the origi- nal least-squares problem (1) with respect to one block of factor variables while keeping all others fixed. Repeating
In
(XC Xn A) A)im An Mi
Ee Wiriz.na ajia'
3

---

Fast Tensor Complction via Approximate Richardson Iteration
this process converges to a nontrivial local optimum and re- duces the original nonconvex problem to series of (convex) linear regression problems in each step (see Appendix A.I).
which is equivalent to solving I, Kronecker regression prob- lems in R/
Updating a block of variables in an ALS step of a TD prob- lem is often a highly structured regression problem that can be solved very fast with problem-specific algorithms. We describe the induced structure of each ALS update for the tensor decomposition types above.
3. Approximate Richardson Iteration
We now present our main techniques for reducing tensor completion to tensor decomposition. When using ALS to solve a TC problem, we must efficiently solve least-squares problems
min Px q
CP factor matrix update = Khatri-Rao regression. In each ALS step for CP decomposition, all factor matrices are fixed except for one, say A (n), ALS solves the following linear cqasresuuares problem:
The rows of the design matrix P correspond to the subset of observations in the TC problem, which means P does not necessarily have the structure of the design matrix in the full TD problem.
N O 1=
A (n) < arg min AER
A
T
x
(6)
A direct approach is to compute the closed-form solution x = (P PP q, but computing (P PP is often im- practical. Two techniques are commonly used to overcome this: (1) iterative methods and (2) row sampling. Iterative methods repeat the same relatively cheap per-step computa- tion many times to approximate the original expensive com- putation. Row sampling methods (e.g.. leverage score sam- pling) randomly pick rows of P and solve a least-squares problem on the sampled rows to obtain an approximate solu- tion to the original problem with high probability, Directly computing leverage scores for general P. however, is also prohibitively expensive since it requires computing the same matrix (P P)-1 (see Appendix A.2 for details).
(m)
Then, we set X, = ||a!5)) 112 for each r E [R] and normalize
the columns of A("), Each row of A (n) can be optimized independently, so (6) solves Jn linear regression problems where the design matrix is a Khatri-Rao product.
Tucker core update- Kronecker regression. The core- tensor ALS update solves the following: for (n) fixed,
N
9 -
aig min
X
A
(n)
vec(9") - vec(
S'ERRI
where the design matrix is a Kronecker product of the fac- tors.
We show that our catroxximaaemmiii--ae method is a prin- cipled approach for tensor completion. In Section 3.1, we prove that lifting restores the structure of the full TD ALS update step, enabling fast least-squares methods for a larger but equivalent problem. In Section 3.2, we show that iter- atively solving the lifted problem (i.e., mini-ALS) is con- nected to an iterative method called the Richardson iter- ation (Richardson, 1911), which we can also view as a matrix-splitiing method. In other words, mini-ALS and the Richardson iteration with a certain preconditioner give the same sequence of iterates (x(*) Jx20- Lastly in Section 3.3, we prove novel convergence guarantees for approximately solving the lifted problem (i.e., for
Tucker factor update = Kronecker-matrix regression. When ALS updates A (n) with all the other factor matrices and core tensor fixed, it solves:
A (x)
arg min AER X Ri,
A
G AT
X
11
where G(n). X(n) are the mode-n unfoldings of 9 and X. This is equivalent to solving J, independent linear regres- sion problems, where the design matrix is the product of a Kronecker product and another matrix. It can be viewed as solving ln instances of structured but constrained linear regression (Fahrbach al., 2022).
mittinattte--miiAii- This allows us to directly use fast leverage sccree sampling al- gorithms for CP decomposition (Cheng et al., 2016; Larsen & Kolda, 2022; Bharadwaj et al., 2023), Tucker decompo- sition (Diao et al., 2019; Fahrbach et al., 2022), and TT decomposition (Bharadwaj et al., 2024) as blackbox subrou- tines. All missing proofs are deferred to Appendix B.
TT-core update = Kronecker regression. Given TT de- and n E [ the left chain Acn E composition (AA (') >
R(h-- R,, and right chain A nn E R x (I. are matrices that depend on the cores A (x') for n' and n' n. respectively (see Appendix A.3 for details). When ALS updates A (i) with all other TT-cores fixed, it solves
3.1. Lifting to a Structured Problem
Consider the linear regression problem with P E R| and q E R given by
A (n)
arg min BER S.pRR-1x110 XR.
(A A TT:
XTT) ||p
x arg min
Px - 2x-l112-
(7)
XER

---

Fast Tensor Complction via Approximate Richardson Iteration
If there exists a tall structured matrix A RIxR with a sub- set of rows 92 [7] such that An P (permutations of the rows allowed), then can lift (7) to higher-dimensional problem while preserving the optimal solution.
Lemma 3.5. Let A, A.P. E RIXR, q E IR! such that P A
and [P q are orthogonal, i.e., Then, the iterative method
(P A)
P =
q = =q (A x
Lemma 3.1. Let b E IR' be the lifted response such that ba = q and byz is a free variable. If
x (1+11) arg min |Ax
1
XERR
(x"
arg min DEjER
||Ax 1 bll2
(8)
simulates Richardson iterations with preconditioner A A for the regression problem minx PPx
912- i.e.
then x" also minimizes (7), i.e., the original linear regres- sion problem min PPx 9l12.
x (++ =x (A A)
P Px (k)
P
(9)
Remark 3.6. In the tensor completion setting, A P van- ishes over 22 so q()) only updates entries in 52 while main- taining q on Ω Thus, computing x ((+11 corresponds to
Lemma 3.2. Problem 8 a convex quadratic program.
Remark 3.3. Problem 8 not a linear regression problem with (structured) design matrix A because there are by vari- ables in the response. However, there is enough structure to employ block minimization to alternate between minimizing x and by-
x (k+1) arg min || - (
XER
arg min {|[Px q12 + +ATT (x- 11
XERR
3.2. Iterative Methods for the Lifted Problem
3.3. Approximately Solving the Lifted Problem
Iterative methods for solving linear systems and regression problems have a long history and have been used to speed up several algorithms in theory and practice. The algorithms we consider use the exact arithmetic model, but all of these methods can be carried out with numbers with log s/s bits, where is the condition number of the matrix (see, e.g., Ghadiri et al. (2023b; 2024)). There is a literature on inexact Richardson iteration for solving linear systems, but they require the error ē to be smaller than than 1/k, which is not achievable with leverage-score sampling (Golub & Overton, 1988; Golub van der Vorst, 1997).
We have shown that alternating minimization for the lifted problem (8) is closely connected to preconditioned Richard- son iteration and inherits its convergence guarantees. For this to be useful, we need to use fast regression algorithms for the (*+11) updates that exploit the structure of A, i.e., when solving minx || - q( [2, where x (k+1) is the vec- tor produced in iteration h of Algorithm 1.
This is where leverage score sampling comes in to play. We exploit the structure of A to efficiently compute its leverage scores, and then we solve the regression problem efficiently but approximately.
Lemma 3.4 (Preconditioned Richardson iteration, (Leee & Vempala, 2024, Lemma 6.))). Consider the least-squares problem x" arg min |Px gll- Let M be matrix such that P M & B P consider the Richardson iteration:
Our next result shows how using approximate least-squares solutions in each step of block minimization affects the con- vergence guarantee of our lifted iterative method. Theorem 3.7. that P and [P q are orthogonal, and
P for some 3 2 1, and
Let A, P
P ERI q E IR,, and B 2 1 such
((++1)= (k) M (x 1M:0--P7 pT
P P'PSAAASB.PT
Then, we have that
Let E $ E [0, 1/33) and approx be an algorithm that for any & E R and q (A-P) computes x E IR in time O(T) such that
[x(++]]- - xx 1=I S (1- 1 [|x x"
- Il.
For the rest of this section, let P R q R' be the zero-masked lifted matrix and vector such that
|AX - f12 < (1 (1+ min Ax - Axl
Then, Algorithm I returns an approximate solution ž E R using approx-least-square squares as a subroutine, such that
Po pT) (As, 0) and (ğa, t = (ba, 0).
2
We now present a key lemma showing that alternating mini- mization between x and by corresponds to preconditioned Richardson iterations on the original least-squares problem. Below, one can easily check that A, P. and q in our lifted approach satisfy this condition.
PP ف| <
1+
min Pxx
q 22
(2/0 Vē
+ PP
PTq
P
2:
in 0
B [-VEB
T Tlog 8/6 time.
5

---

Fast Tensor Complction via Approximate Richardson Iteration
Algorithm 1: approx-mini-als
algorithm for approximately solving the structured least- squares problem in each step of approximate miniatee-mmiiiSS... For given E E (0, 1/89), the approximate solver computes solution X E R in time O(Ta) such that
Data: A, A, 1, E (0,1), E 1/ with P A 33.PTP Result: X EE
Ax b S (1+ min Ax bll..
Initialize (i)iizexx())_
for =0, m0,,==0,,, 0,I,....; (2% do
21//-
Therefore, for desired e E (0, 1), set ( and use a sufficiently small & t E2 in Theorem 3.7. Putting everything together, Algorithm finds an approximate solu- tion E IR in time O 518 llg 8 that satisfies
Set q()) + q (A (*)
// implicit
Set x(++1) to vector such that
|Ax S (1+0) ming | -x--xx0 1
52
return (A)
Px - q11 S (1+ I| qll + Teql:: (10)
where Tp and Tp the orthogonal projection matrices into the column space and null space of P, respectively (see Appendix A. With this in hand, we are ready to present the running times of our lifted iterative method for TC prob- lems by combining Theorem 3.7 with state-of-the of-the-art tensor decomposition results based on leverage score sampling.
Remark 3.8. To better understand Theorem 3.7, observe that P P P P is a 9-spectral approximation of A A, & controls the reducible error & |Px I|3. and (1 + E) is the error in the approximate least-square update for each (k))
Bounding 3. First, observe that in the case of TD, we have = A, so 3 1. More generally, if rank () = x min R} and A UΣVT is a compressed SVD, then A is said to satisfy the standard incoherence condition with parameter μs (Chen, 2015)
4.1. CP Completion
Each ALS update step for CP completion solves a regression problem where the design matrix the Khatri-Rao product:
for A (k) IR R
A #k
O
N A (n) E g
max sE
e TU S
uss I
max PEIR
V e,
1l
us 5
and Q = (X (k) E RR
R
The l|e_ U12 and T e,.11 values are the leverage scores of the rows and columns of A. Applying Cohen al. (2015, Lemma 4), if each row of A is observed with probability p que logs for some absolute constant c, then
A (k) + arg min
(A In A QQ/E
AERAXR
The design matrix (A *k) 12 does not necessarily have any structure, so direct method relies on solving the normal equation, which takes O R[ + Ix)) time. Thus, the running time of one round of ALS, i.e., updating all N factors, is O(N(R* + R (1) + R( Σ Tn).
such that 2
3 A ATA ATA,
AT ATA
which gives 3 = 2/p. Let Ç maxie(] |ajll2, where a, is row i of A. Then, the at² rid-redge leverage scores of A (i.e., a (A a Ix)) - for 2 are at most 1/a. If is the observation rate, taking = clogs gives the required incoherence condition. This can be done by introducing an fz-regularization term to the TC optimization problem (i.e., solving a ridge regression problem each ALS step). Note that usually a can be chosen to be much smaller in practice.
Previous work on CP tensor decomposition (Cheng et al., 2016; Larsen Kolda, 2022; Bharadwaj et al., 2023) de- veloped fast methods for efficiently computing the leverage scores of Khatri- product matrix. particular, Bharad- waj al. (2023) designed data structure for computing and maintaining the leverage scores of A # during ALS updates. This approach requires sampling Ö(R/s) rows of A A** Due to the Khatri-Rao product structure, each row of A # can be mapped to sequence of one choice from the rows of each (n) for n E [N]](k) Hence, sam- pling row from A* is equivalent the following: for each n [[]] sample a row from A (n)) according to some conditional distribution given sampled rows from A (1) .... A (n-1), and then compute the Hadamard product of N sampled rows. Maintaining the full L, -i-dimensional vector for a conditional probability for each n costly, so Bharadwaj et al. (2023) developed binary tree- bre-based data structure to speed up leverage-scon score sampling for A #
4. Sampling Methods for Tensor Completion
We are now ready to efficiently solve the unstructured least- squares problem (7) induced by ALS for tensor completion, ie., for P E IR R| R and observations q E [R find
x" arg min PXx all2.
X
As in Algorithm 1, we lift this problem to higher dimension to get structured design matrix A, and use a known fast
6

---

Fast Tensor Complction via Approximate Richardson Iteration
Applying Bharadwaj et al. (2023, Corollary 3.3), round of ALS runs in time Ō(s-1 Ex= (T, R + RNR))) Using their CP TD algorithm the approximate solver in Algo- rithm 1, and combining its guarantee with Theorem 3.7, we can extend their approach to CP completion.
4.2.2. FACTOR MATRIX UPDATE
The ALS factor matrix update for A (k) in the Tucker com- pletion problem is
N (nn)
Corollary 4.1. There is an ALS CP completion algorithm such that (i) after a factor matrix update, each row of A ( satisfies (10), and (ii) the total running time of one round is
A (K) f arg min AER K
(n)
G (k)
A T Q
11 n.lk
$2
where Q = (X (k) E RIjKK is sparse matrix of ob-
3 of
N
1
servations. The running time of a direct method that solves the normal equation is ORR + Rx| + Rs + Ik)), where Ryx = Rx-
Σ E1 R=1
RT
NR log
Note that there is no dependence on [S2] in the running time due to leverage score sampling, i.e., it runs in sublinear time.
The running time of the sampling- factor-matrix up- date for A (R) in Fahrbach al. (2022, Theorem 1.2) for the full decomposition problem is
4.2. Tucker Completion
N
N e2
+ R2
Fahrbach et al. (2022) designed block-sketching techniques and fast Kronecker product-matris multiplication algorithms to exploit the ALS structure for Tucker decomposition.
Rn +
+ RRRR
Combining this result with Theorem 3.7, Algorithm has the following running time for factor matrix update.
4.2.1. 1. CORE TENSOR UPDATE
Corollary 4.3. There is an algorithm that computes an ALS
Recall that for Tucker decomposition we use the notation I [N In and R [I Ru- The ALS core tensor update in the Tucker completion problem is
Tucker completion factor matrix update for A (k), with each row of A (k) satisfying (10), in time
o o((ct
82-
+ 1RE n== Rn Blog
A (n)
9 -
arg min
veo(S') - vec(X)
61
S'ERR,X XRN
n=11
Ω
where C Em Ru + g4 N2
(In
The design matrix above restricted to Ω is exactly P in our general setup.
4.3. TT Completion
We compare the running times of the direct method and our lifting approach. In the former, we can compute an exact solution to least-squares problem x* = P) in time ( R2 1 ).).
Each ALS step for TT decomposition solves the following seastesquares problem with a Kronecker product-type design
matrix: for A# == Ack & A R. lu x (R, R and Q = (k)9 E IR Rix
To achieve fast lifted method, we solve the second step of Algorithm using the leverage score sampling-based core tensor update algorithm in (Fahrbach et al., 2022, Theorem 1.2) with running time
A (k) arg min
(A (B QlIF
BE Re-10- X
Solving this directly with the normal equation takes O + R& Ru + T&)) time for R3 == R RA. Thus, the time for round of ALS is O - [Rg R2 ( (R,, In))).
N2 E²
R2-0"
Σ Rn
E
In contrast, Bharadwaj et al. (2024, Corollary 4.4) show that one round of approximate TT-core updates, if R, = R for all n E [N -1]. can run in time Ö(R*s Σm (N + I.)). See Appendix C for technical details.
where 0* is an optimizable constant depending on [Ra)eEJN] INE[N-- Using this the opprox-last--quarre sub- routine in Theorem 3.7, achieve the following.
Corollary 4.2. There is an algorithm that computes an ALS Tucker completion core tensor update satisfying (10) in time
Corollary 4.4. There is an ALS TT completion algorithm
such that (i) after a TT-core update, each row fiber of (m)) satisfies (10), and (ii) the total running time of one round
((c+
88 8-22-00
1
3 log
(11)
E1
E2
N (N + In) log
δ
Σ n=1
1
where C:= EN= (I + +.8 N
82
7

---

Fast Tensor Complction via Approximate Richardson Iteration
5. Experiments
In this section, we study the empirical performance of our algorithm and compare it to direct and expectation maxi- mization (EM) methods for a coupled matrix problem and CP completion tasks.'
5.1. Warm-Up: Coupled Matrix Problem
First we consider the coupled matrix problem
Figure I. Coupled matrix results for E R=**, X, Y E R Rdxd with 2000 and d 10 that compare the direct method, mini- via leverage score sampling.
AXB + CYD = E,
ALS, and approximate-minie
where A, L,B,C.D C. E Raxd are given, X, Y E Rdddd are unknown, and E E RRxx is a matrix with half of its entries randomly revealed (Baksalary & Kala, 1980). For fixed Y, we compute X by solving the Kronecker regression problem
does not necessarily attain the minimal error in the matrix's s kernel space. i.e., the first term on the right-hand side of (10) can exceed the minimum error by a factor of (1 + e), which allows it to follow a different convergence path and achieve a lower MSE than the other methods around step 15 of ALS.
arg min| | (Bog % A) vee(X) vec(E CYD Il22. (12) XERaxa
The matrix Y can be updated by solving similar regression problem. Therefore, we can apply an alternating minimiza- tion algorithm to compute X and Y. For these experiments, we initialize X Y = I. We present the results in Figure I, which are averaged over five trials.
5.2. CP Completion
Now we compare methods for the CP completion task. For given tensor X and sample ratio p E [0, 1], let Xg be a partially observed tensor with a random p fraction of entries revealed. We fit Xg with a rank-R CP decomposition by minimizing the training relative reconstruction error (RRE) | ((C-X)allee x/ [|Caele using different ALS algorithms.
Data generation. The entries of X, Y are sampled independently from uniform distribution on [0, 1), and then we set E + AXB T + CYDT. We consider the setting where half of the entries of E are observed (chosen uniformly at random). We set n 2000 and d = 10. Note that since we observe a subset of entries of E, the Kronecker regression structure lost.
Datasets. We consider two real-world tensors. CARDIAC- MRI is a 256 x 256 x4 tensor of MRI measurements in- dexed by (a., y, 2, t) where (z, y. 2) is a point in space and t corresponds to time, HYPERSPECTRAL is 1024 x 1344 x 33 tensor of time-lapse hyperspectral radiance images (Nasci- mento et al., 2016). We also consider synthetic low-rank CP and Tucker tensors in Appendix D.1.1.
Algorithms. We compare the direct, mini-als, and approximate-mini- methods. The latter two use adap- tive step sizes based on the trajectory of the iterates (*) (see Appendix 1.1.2 2.222 for details). direct solves the normal equa- tion in each ALS step and runs in (( +R²) time since mini-als is Algorithm 1 with = and e > 0. mini mini-als uses the Kronecker product proper- (B B (A (A)A) and vec vec(AXB) (AXB for improved efficiency. approximater -iiimmiiiiii nals uses leverage score sampling for Kronecker products similar to Fahrbach et al. (2022); Diao et al. (2019), which is a direct application of Corollary 4.2. For leverage score sampling, we sample 1% of rows in each iteration of wwwiatt-maaiiiiiieeeeeeeiee
We compare direct, parafac, mini-als, and accelerated -iii-minn methods. direct solves the normal equation in each ALS step for the original problem of the form (4) and has rumning time ( + R³) since it computes (A An) mini-als is our lifting approach in Algorithm with ē = O and Æ > O. parafac is the EM algorithm of Tomasi & Bro (2005), which is equivalent to running Algorithm I for exactly one iteration each step of the ALS algorithm. accelerated iiieeed--ii--ii -als uses adaptive step sizes based on the trajectory of the iterates (8)
lgorittsms.
computes (P PP)
ties (BB 8 A) (B @ A))
(B @ A) vec(X
Results. The left plot in Figure shows the total running time of these three algorithms across all ALS iterations. The right plots shows the mean squared error (MSE) at each step of ALS. Since approximate-miini-ll is stochastic, it
Results. In Figure 2, we present plots for CARDIAC- in the top row and HYPERSPECTRAL. in the bottom row.
In the first column, we set p = 0.1 and sweep over the rank R using di direct demonstrate how the train RRE (solid line) and test RRE (080X//r//ll ||]0/// (dashed line)
Thhe code is available https://githhee com htt.uu...comaaabaach// fast- caansor-comletioon
8

---

Fast Tensor Complction via Approximate Richardson Iteration
21
2400
250
144
7000
1.00
Figure Algorithm comparison for a low-rank CP completion task on CARDIAC-MRI HYPERSPECTRAL tensor datasets. The first column illustrates convergence rates and relative reconstruction error (RRE) of the ALS algorithm for various CP-decomposition ranks, The second column plots RRE over ALS steps for direct, parafac, and our method (mini-als) under different values. The third column shows total running time of these algorithms varying sample ratios, i.c., number of observed entries. Finally, the fourth column displays total running times for direct, parafac, and accelerated -iin-als (instead mini-als) across different & settings.
decrease as function of R and the ALS step count. Note that the test RRE the loss on the entire (mostly unobserved) tensor, and slightly above the train RRE curves all times.
illustrate this idea with Figure in Appendix D.1.2. accelerated-mini achieves better solution qual- ity than mini-als for given and always runs faster, especially for small values of p.
In the second column, we fix the parameters (p, R) (0.1, 16) and study the solution quality of m'mini-als for different values of e compared to di direct and parafac. As e decreases, we recover solutions in each step of ALS that are as good as direct, which agrees with the claim that lifting simulates the preconditioned Richard- son iteration (Lemma 3.5).
6. Conclusion
This work introduces novel lifting approach for tensor com- pletion. We build on fast sketching- TD algorithms as ALS subroutines, extending their guarantees the TC set- ting and establishing novel connections iterative methods. We prove guarantees for the convergence rate of an approxi- mate version of the Richardson iteration, and we study how these algorithms perform in practice on real-world tensors. One interesting future direction to analyze the speedup that adaptive step sizes give accelerated-min -als.
In the third column, we sweep over for R = 16 and plot the total running time direct, parafac, and mini-als for 10 rounds of ALS. The running time of direct increases linearly in p [||] as expected. In- terestingly, the running time of mini-als decreases [S2] increases, for < 0.1, since the lifted (structured) matrix AT A becomes a better preconditioner for the TC design matrix Ag (i.e., 3 + as + 1). This means mini-als needs fewer iterations to converge. Finally, parafac performs one mini-ALS iteration in each ALS step and is therefore faster but converges at slower rate.
A Acknowledgements
We would like to thank Arkadi Nemirovski and Devavrat Shah for insightful discussions. MG and AJ acknowledge the support of ARO MURI W911-NF- 19-N-1991-- 9-1-0217 and ONR N00014-23-1-2299 YK supported in part by NSF Award CCF-2106444, and did this work student researcher Google Research.
In the fourth column, compare the total running algorithm (dash-
time of the accelerated-ming dot lines) to mini-ali (third column). Our accelerated method extrapolates the trajectory of the iterates (() during mini-ALS using a geometric series, which al- lows us solve collinear iterates in one step. We
Impact Statement
This paper presents work whose goal is to advance the field of machine learning. There are many potential consequences of our work, but none that we feel must be highlighted here.

---

Fast Tensor Complction via Approximate Richardson Iteration
References
Fazel, M. Mairix Rank Minimization with Applications.
PhD thesis, Stanford University, 2002.
Acar, E., Dunlavy, D. M., Kolda, T. G., and Mørup, M. Scal- able tensor factorizations for incomplete data. Chemomet- rics and Intelligent Laboratory Systems, 106(1):41-56, 2011.
Filipović, M. and Jukić, A. Tucker factorization with miss-
ing data with application to low-n-rank tensor completion. Multidimensional Systems and Signal Processing, 26(3): 677-692, 2015.
Bader, B. W. and Kolda, T.G. G. Tensor toolbox for MATLAB, version 3.6. https: //www. tensortoolbox. org/, 2023.
Ghadiri, M., Fahrbach, M., Fu, G., and Mirrokni, V. Approx- imately optimal core shapes for tensor decompositions, International Conference on Machine Learning, pp-
Baksalary, J. and Kala, R. The matrix equation A CYD E. Linear Algebra and its Applica- tions, 30:1:00:444 147, 1980.
11237-11254, PMLR, 2023a.
Ghadiri, M., Peng, R., and Vempala, S. S. The bit complex- ity of efficient continuous optimization. 2023 IEEE 64th Annual Symposium on Foundations of Computer Science, pp. 2059-2070, IEEE, 2023b.
Barak, B. and Moitra, A. Noisy tensor completion via the sum-of-squares hierarchy. In Conference Learning Theeory pp. 417-445. PMLR, 2016.
Bharadwaj, V., Malik, O. A., Murray, R., Grigori, L., Buluc, A., and Demmel, J. Fast exact leverage score sampling from Khatri-Rao products with applications to tensor decomposition. In Advances in Neural Information Pro- cessing Systems, volume 36, pp. 47874 1877447991, 2023.
Ghadiri, M., Lee, T.,, Padmanabhan, S., Swartworth, W., Woodruff, D. P., and Ye, G. Improving the bit complexity of communication for distributed convex optimization. In Proceedings of the 56th Anmual ACM Symposium on Theory of Computing, pp. 1130-1140, 2024.
Bharadwaj. V., Rakhshan, T., Malik, O. A., and Rabusseau, G. Efficient leverage score sampling for ten- sor train decomposition. In Advances in Neural Informa- tion Processing Systems, 2024.
Golub, G. and Overton, M. L. The convergence inexact Chebyshev and Richardson iterative methods for solving linear systems. Numerische Mathematik, 53(5):571-5953 1988.
Golub, G. H. and van der Vorst, H. A. Closer to the solution: Iterative linear solvers. Institute of Mathematics and its Applications Conference Series, volume 63, pp. 63-92. Oxford University Press, 1997.
Candes, E. and Recht, B. Exact matrix completion via convex optimization. Communications of the ACM, 55 (6):111-119, 2012.
Chen, Y. Incoherence-optime matrix completion. IEEE Transactions on Information Theory, 61(5):2909-2923, 2015.
Harris, C. R., Millman, K. J., van der Walt, S. Gommers, R., Virtanen, P., Cournapeau, D., Wieser, E., Taylor, J., Berg. S., Smith, N. J., Kern, R., Picus, M., Hoyer, S., van Kerkwijk, M. H., Brett, M., Haldane, A., del Río, J. F, Wiebe, M., Peterson, P., Gérard-Marchang P., Sheppard, K., Reddy, T., Weckesser, W., Abbasi, H., Gohlke, C., and Oliphant, T. E. Array programming with NumPy. Nature, 585(7825) 2020.
Cheng, D., Peng, R., Liu, Y., and Perros, SPALS: Fast alternating least squares via implicit leverage scores pling. Advances in Neural Information Processing Sys- tems, 29, 2016.
Cohen, M. B., Lee, Y. T., Musco, C., Musco, C., Peng, R., and Sidford, A. Uniform sampling for matrix approxi- mation. In Proceedings of the 2015 Conference on Inno- vations Theoretical Computer Science, pp. 181-190, 2015.
Healy, M. and Westmacott, M. Missing values in experi- ments analysed automatic computers. Journal the Royal Statistical Society: Series C (Applied Statistics), 5 (3):203-206, 1956.
Hillar, C. J. and Lim, L.-H. Most tensor problems are NP- hard. Journal of the ACM, 60(6): 2013.
Diao, H., Jayaram, R., Song. Z., Sun, W., and Woodruff, D. Optimal sketching for Kronecker product regression and low rank approximation. Advances in Neural Information Processing Systems, 32, 2019.
Jain, P. and Oh, S. Provable tensor factorization with miss- ing data. Advances in Neural Information Processing Systems, 27, 2014.
Fahrbach, M., Fu, G., and Ghadiri, M. Subquadratic Kro- necker regression with applications to tensor decomposi- tion. Advances in Neural Information Processing Systems, 35:28776-287899 2022.
Kilmer, M. E. and Martin, C. D. Factorization strategies for third-order tensors. Linear Algebra and its Applications, 435(3):641 1658, 2011.
10

---

Fast Tensor Complction via Approximate Richardson Iteration
Kolda, T. G. and Bader, B. W. Tensor decompositions and applications. SIAM Review, 51(3): 2009.
Song, Q., Ge, H., Caverlee, J., and Hu, X. Tensor comple- tion algorithms in big data analytics. ACM Transactions Knowledge Discovery from Data, 13(1):1-48, 2019.
Kossaifi, J., Panagakis, Y., Anandkumar, A., and Pantic, M.
TensorLy: Tensor learning in Python. Journal of Machine Learning Research (JMLR). 20(26), 2019.
Tarzanagh, D.
and Michailidis, G. Fast randomized algorithms for t-product based tensor operations and de- compositions with applications imaging data. SIAM Journal Imaging Sciences, 11(4):2629-2644, 2018.
Larsen, B. W. and Kolda, G. Practical leverage-based sam- pling for low-rank tensor decomposition. SIAM Journal on Matrix Analysis and Applications, 43(3):1488-1517, 2022.
Tomasi, G. and Bro, R. PARAFAC and missing values. Chemometrics and Intelligent Laboratory Systems, 75(2): 163-180, 2005.
Lee, Y. and Vempala, S. Techniques in Optimization and Sampling. 2024. URL https: //github. com/YinTat/ optimizationbocke http///boooo//m//iiim/nnaa// pdf.
Wang. M., Yu, Y., and Li, H. Randomized tensor wheel decomposition. SIAM Journal on Scientific Computing, 46(3):A1714-A1744 2024.
Little, R. J. and Rubin, D. B. Statistical Analysis with Missing Data, volume 793. John Wiley Sons, 2019.
Wu, Z.-C., Huang, T.Z.., Deng, L.J... Dou, H.-X., and Meng, D. Tensor wheel decomposition and its tensor completion application. Advances in Neural Information Processing Systems, 35:27008-27063 2022.
Liu, A. and Moitra, A. Tensor completion made practical. Advances in Neural Information Processing Systems, 33: 18905-18916, 2020.
Zhao, Q., Zhou, G., Xie. S., Zhang, L., and Cichocki, Tensor ring decomposition. arXiv preprint urXiv:1606.05533 2016.
Malik, O. A. and Becker, S. A sampling-based based method for tensor ring decomposition. In International Conference on Machine Learning, pp. 7400-7411, PMLR, 2021.
Zheng, Y.-B., Huang, T.-Z., Zhao, X.-L.,, Zhao Q., and Jiang, T.-X. Fully-connected tensor network decomposi- tion and its application to higher-order tensor completion. Proceedings of the AAAI Conference on Artificial In- telligence, volume 35, pp. 11071-11078 2021.
Malik, O. A., Bharadwaj. V., and Murray, R. Sampling- based decomposition algorithms for arbitrary tensor net- works. arXiv preprint arXiv:22/0.03822 2022.
Montanari, A. and Sun, N. Spectral algorithms for ten- sor completion. Communications on Pure and Applied Mathematics, 71 71(11):2381-2425, 2018.
Nascimento, S. M., Amano, K., and Foster, D. H. Spatial distributions of local illumination color in natural scenes. Vision Research, 120:39-44 2016.
Oseledets, I. V. Tensor-train decomposition. SIAM Journal Scientific Computing, 33(5):2255-23177 201
Richardson, L. F. The approximate arithmetical solution by finite differences of physical problems involving dif- ferential equations, with an application the stresses in a masonry dam. Philosophical Transactions of the Royal Society of London. Series A, Containing Papers of a Mathematical or Physical Character, 210(459-470): 307-357, 1911.
Shah, D. and Yu, L. Iterative collaborative filtering for sparse noisy tensor estimation. In 2019 IEEE Inter- national Symposium on Information Theory, pp. 41-45. IEEE, 2019.
Shah, D. and Yu, C. L. Robust max entrywise eπror bounds for tensor estimation from sparse observations via similarity based collaborative filtering. IEEE Trans- actions on Information Theory, 6955)3121-33 2023.

---

Fast Tensor Complction via Approximate Richardson Iteration
A. Missing Details for Section 2
A.1. Least-Squares Linear Regression
For a design matrix A E R2xd and response b E IR", consider the least-squares problem
x" arg min Ax - bll2.
XER
The solution x" can be obtained by solving the normal equation A T Axx A A Therefore, x (ATA) AA A A is singular, then we use the pseudoinverse A
The orthogonal projection matrix : RR" + R onto the image space of A is defined by TA = (ATA))-AAT and satisfies A = and || < livi2 for any E R.. Recall that any E R" can be uniquely decomposed as V = = + AA V, where TA = I, T is the orthogonal projection to the orthogonal subspace colsp(A).
Given AA b + TA i b, the first term is the reducible error by regressing b on x, i.e., taking the optimum x so that Ax" T.b. b. The second term A is the irreducible error, i.c., min || bll2 = A||2.
A.2. Leverage Score Sampling for Tensor Decomposition
ALS formulations show how each tensor decomposition step reduces to solving a least-squares problem the form minx || bll2 with highly structured A. While we can find the optimum in closed form via (A A) A T matrix A has J = ..... Ty rows corresponding each entry of the tensor (i.e., it is tall skinny matrix), which can make using the normal equation challenging in practice.
Randomized sketching methods are popular approach to approximately solving this problem with faster running times with high probability. In general, these approach sample rows of A according to the probability distribution defined by the leverage scores of rows, resulting in random sketching matrix S whose height much smaller than that A. For a matrix >> R), the leverage scores of A the vector E [0, 1}' defined by
A E R with
4 def (A AT) AT....
Then, for given e, 8, E (0,1), the sketching algorithm samples ((*/26) many rows, where the i-th row drawn with guarantee that
probability 4 Σ, f./rank (A). With probability least 1 5,
min SAA Sbll2 s (1 + e) min Ax - bll2-
The reduced number of rows in SA leads to better rumning times for the least-squares solves. However, naïvely computing leverage scores takes as long computing the closed-form optimum since need to compute (A AA)t.. This is where the structure of the design matrix A comes in to play, i.e., to speed up the leverage score computations.
A.3. Tensor- Train Decomposition
Given a tensor-train (TT) decomposition {A (A( Jn=1 and index E [N], define the left chain Acn (1)xx= and the right chain A E R.R. (- as:
Gen .... 1n-1-
Σ
II
Po k=1
N
azn("nn Antt
Σ
II
Ta-
Fwo k=n+1
where for any is E [44] with s(+ n) [N], i i--1 == 1+ 11 and int1 in == 1
1) j=nt1 I,. When ALS optimizes A (x) with all other TT-cores fixed, it solves the regression problem:
A (m) arg min
(A @ AAJ. BT) X 1n
BE xx
which is equivalent to solving I, Kronecker regression problems in Th. Ill
12

---

Fast Tensor Complction via Approximate Richardson Iteration
A.4. Tensor Networks
tensor network (TN) is powerful framework that can represent any factorization of tensor, so it can recover the three decompositions above as special cases. A TN decomposition TN(A ................................ A (N)) represents given tensor X with N tensors (1) ;.....A' (M)
A
A(M) and a tensor diagram. As in the above decompositions, the goal is to compute
,...,, arg min - TN(A ....., A
tensor diagram" consists of nodes with dangling edges, where a node indicates a tensor, and its dangling edge represents mode, so that the number of dangling edges is the order of the tensor. For example, a node without an edge indicates a scalar, one with one dangling edge is vector, and one with two dangling edges is matrix.
When two dangling edges of two nodes are connected, we say that the two tensors are contracted along that mode (i.e., mode product of those two tensors). For example, when a node with two dangling edges shares one edge with another node with one dangling edge, it indicates a matrix- vector multiplication. Hence, the number of unmatched dangling edges in a tensor diagram corresponds to the order of its representing tensor.
notable special case a fully-connected tensor network (FCTN) (Zheng et al., 2021) decomposition. It consists of N tensors of the same order, where in its tensor diagram, all pairs of nodes are connected as its name suggests.
ALS for TN decomposition. Given TN decomposition {A((n) (A(()))E(N) fixed, it solves linear regression problem of the form
JNE[M), when ALS optimizes tensor A (n) with all others
A arg min A B- X||r,
B
where Ayn is an appropriate matricization depending on A A (n-1), ( (x+11 matricizations of B and X, respectively. Structure of A can be specified through a new tensor diagram obtained by removing the node of A (n)) from the original tensor network diagram.
(N), and B and X are suitable
Just as the ALS approaches for other decomposition algorithms, Malik et al. (2022) proposed a sampling- based approach via leverage scores. First of all, they pointed out that A *n A can be efficiently computed by exploiting inherent structure of A (i.c., contract a series of matched edges a tensor diagram in an appropriate order). They then presented a leverage-score sampling method that draws rows of y without materializing a full probability vector, and in spirit this approach is similar to one used for the CP decomposition in Section 4.1.
Other tensor decompositions. We also discuss other important special cases of TN decompositions including the tensor ring (TR) (Zhao al., 2016) and tensor wheel (TW) (Wu al., 2022) decompositions. Both decompositions can be succinctly described through tensor diagram notation.
where for each & E [N], the node of (A) is connected to the two neighboring nodes of A (k1) and A (k+1) (precisely, the superscripts here refer to modulo with respect to N) that the resulting tensor diagram looks exactly like ring. Since fast TR decomposition through leverage-soon score sampling is facilitated by Malik & Becker (2021), our approach of lifting and completion can be straightforwardly applied to the TR decomposition.
A 'TR decomposition consists of N tensors ( ((N), ( (N),
(k)
TW decomposition consists of one core tensor and N factor tensors (1) A (N), As its name suggests, its tensor diagram is essentially that of the tensor ring (consisting of N factor tensors) with the node representing e connected to all factor tensors. Thanks to recent work of Wang et al. (2024) using leverage-score sampling for the TW decomposition, approach can be readily applied as well.
The t-SVD (Kilmer& & Martin, 201 1) and I-CUR (Tarzanagh & Michailidis, 2018) decompositions for third-order tensors that have been applied to spatiotemporal data. They factor a tensor via the t-product into three tensors; in the case of t-SVD, one of these factors is f-diagonal (see Kilmer & Martin (2011) for definitions). Tarzanagh & Michailidis (2018) propose fast, sampling-bast based algorithms for t-product-
tensor decompositions. These algorithms can be combined with our lifting method and thus extended to the tensor completion setting.
We refer readers to http:://tennsrreewoorrrr///////e
for more details.
13

---

Fast Tensor Complction via Approximate Richardson Iteration
B. Missing Analysis for Section 3
B.1. Proof Lemma 3.1
Lemma 3.1. Let E R! be the lifted response such that by = q and byy is a free variable. If
(x",
arg min xER by ERI-10
Ax - b
(8)
then x* also minimizes (7), i.e., the original linear regression problem minxERR Px 9/12.
Proof. For any x, have
-Ax- 1 bl12 || bullă + || b
so
min Px
x.b q112 S
min Ax - bll2
Moreover, for any x, taking by
Apx gives |
b
1
0, which implies that
min Px
q
min
Ax
bll2.
x x,b
Therefore,
min Px
q1b
min x,b
mni||Ax Ax
bll2.
and x" also minimizes (7).
B.2. Proof of Lemma 3.2
Lemma 3.2. Problem 8 is convex quadratic program.
Proof. Since q E R defined as qn by and qा
O, we can write (8) the following equivalent manner:
(x", b
arg min XER" b
A
x
b
2
where I is the I identity matrix.
B.3. Proof of Lemma 3.5
Lemma 3.5. Let A, A, RIRR q R such that P A and P are orthogonal, i.e., (P-A)) T ql 0. Then, the iterative method
(()) - (A P) (k)
x (++1) arg min Ax q 11
xERR
simulates Richardson iterations with preconditioner A T A for the regression problem ming PP - -Q13. i.e.
(1+ =x-11 xx() (A AA P -
(9)
Proof. Assume that A ATA A is full rank. Solving the normal equation for x (N+1).
x (A++1)
(A A) -1 A q
= =(A TATA) A) -1 (q + (A P)
(A T - A (Px -
14

---

Fast Tensor Complction via Approximate Richardson Iteration
Since P A and [P
q are orthogonal,
A (Px (8) - ① = (P- (P A)) [P qं xx
--1
= (Px q
Therefore,
x - Px p
which completes the proof.
B.4. Proof of Theorem 3.7
Theorem 3.7. Let A, P E RIxR, qE GERY and 3 2 1 such that P A and P q are orthogonal, and
P PTP.
LetE E (0,1), E (0.1/3²) and squares be an algorithm that for any & E R and f = q (A -P)x (A-P)&
computes X E R in time ((T) such that
|A
fl2 S (1 + 8) min |AAx - ff12.
Then, Algorithm/ / returns approximate solution X R. using seprox-least-squee
subroutine, such that
2
|PX P -
1
min Px 4|p
(1/8
vēy²
+ 5
P
P
1 9
in 0 01-9/0 log 8/6) time.
Proof. We show that Algorithm gives the desired output. Suppose approx-least-squares yields x (l+1) for given inputs A, P, q. and (()) (i.e., & +- (1) &- (*)), and x +- x(*+1)) which satisfies
|Ax -6 1| < (1 + 3) min| Ax q
(()) 1 (1 2) || - 13
We can also decompose the LHS using k(k) = T (k) + + q as follows:
(x+1+1) - ((t) = ||Ax
-
TAA 12 + q
Combining the above, we get
AA (++1) π (k) 1 ~ 13.
Denoting x P P q arg min
Bx - -qll2 and using the triangle inequality,
[x+++ Ax"
< -
- ||a + q Ax
< < (k) 112 + ( Ax
(13)
We now bound each term in the RHS. As for the second term, since q q + (A x))()) and (P (P-A)T
P
q
by Lemma 3.5,
A A) ATQ arg min || q(k) 12
xx() (A (ATA)
Px TT
15

---

Fast Tensor Complction via Approximate Richardson Iteration
which is exactly Richardson iteration with preconditioner M <- A A in Lemma 3.4 (satisfying P ATA
PPP). Thus, [I(A ((AAA))-ATTO x" A
-1-8-11 ||x(() x x"IlAT A and
T Ax 112 S 1 (1-3) }) ||xx(*) Ax" 112
(14)
Regarding the first term in (13), since (x
( is in the column space of A,
TA- (q + (A B) -A- (q- (P)(
T
Therefore,
|| 1 - l Px
= |P - Px (k) 13 + min PP - q
$ Axx A +
Px q
X
where the last inequality follows from PTP 3 ATA. A. Thus,
T ( |Ax (x min PPx q 2.
(15)
Combining (13), (14), and (15), we have
|xxx++1 (k+1) Ax S 1 (1- g vē |AAx Ax + VE min || 71l2.
Denoting a 1 + Vε, by induction, we have
|Ax( Ax" 1/2 S a ||Ax* - Ax 1|2 + (1+ a a + ................................ +a +.... x v min Px - E/ا2
a
at ||xx - Ax + 1-a a
Axx
x VB min Px লlle. 2.
(16)
We also have
PP - -Ql3 = P T
2
PT q l
= Px (k)
Px
min
Px
q
(17)
We then bound the first term by using P 3 A 3 P P and (16) follows:
( - Px
11 < |x (k) Ax" 1
1 # 2
< 2a 2k
Axx - Ax (0) 12+2
E min Px - âl12
1 -
a
2
S 2,3 2h
|Xx Px 11 +
x
min
Px E|p
1 0
x
Putting this bound back into (17),
[x (k)
9 < 230 2k Px (0 +
( (1+28 28
1 - a*
min
Px -q13.
1
Setting
log(23/3
k:
22/33 v²)
we have
12 ( 1+ (1/3 - V2y²
2
P - - <
|Px
(x
min Px 9 2
which completes the proof with (0) 0.
16

---

Fast Tensor Complction via Approximate Richardson Iteration
C. Additional Details for Section 4
Bharadwaj et al. (2024) proposed sampling-basee ALS algorithm that relies a canonical form the decomposition with respect the index k. Any TT decomposition can be converted to this form through a QR decomposition, which ensures that (AT A# IRR-1RA Rg. It follows that the leverage scores A #k are simply the diagonal entries of
th (A (Ac 2Rd @ (A AAA)
follows from properties of the Kronecker product that
Li (A #k) lick (Ack) lizs (A
where j i.. 1 iN. ick i, ... 1- and i> 1..1--12 iN (see Appendix A.3 for definition). Therefore. efficient leverage score sampling for A # reduces to that for Ack and A₃k- To this end, Bharadwaj et al. (2024) adopt an approach similar to Bharadwaj et al. (2023) for leverage score-based CP decomposition. Each row Ac* corresponds to a series of one slice for each third-order tensor A (n) for < k, which results in a series of conditional sampling steps using a data structure adapted from the one used for CP decomposition.
D. A dditional Details for Section 5
All experiments are implemented with NumPy (Harris et al., 2020) and Tensorly (Kossaifi et al., 2019) on an Apple M2 chip with GB of RAM.
D.1. CP Completion
D.1.1. SYNTHETIC TENSORS
We run the same set of experiments as in Section 5 on two different random low-rank tensors:
RANDOM-CP is 100 100 100 tensor formed by reconstructing random rank-16 CP decomposition.
RANDOM-TUCKER is 100 100 tensor formed by reconstructing a random rank-(4, (4, 4, 4) Tucker decomposition.
randor AUCKET
rendor SUCKeT
400 125
1/01
Figure 3. Algorithm comparison for low-rank CP completion task the RANDOM-CP and RANDOM-TUCKE: tensor datasets.
D.1.2. ACCELERATED METHODS
We explain how to accelerate the Richardson iteration.
17

---

Fast Tensor Complction via Approximate Richardson Iteration
(a(1)
(a
200
(a 01
17
(x
(x
Figure 4. In the one marviaiable case, our approach proceeds the following order: (x (/)) + (xx ((0)) - (1) (
1. For odd iterations (e.g., the first iteration), run mini-ALS normally.
2. For even iterations, compute (k++1) using normal mini-ALS, but then set
x(++1) (k) 1 a
1 (((+1) -x(1))-((()) x
-x(0 2 (k))
[
where
Note that setting = O is equivalent to running mini-ALS normally. We explain the intuition behind this accelration with the example in Figure 4, which illustrates the case where æ and b/ both have only one variable. As mentioned previously, the lifted problem is a convex quadratic problem (see Lemma 3.2), and the iterations of mini-ALS alternate between optimizing x and optimizing byg. In Figure 4, the star (*) denotes the optimal point, and the ellipses denote the level sets of the quadratic function. In the steps where we optimize x, we search for the point on the line that parallel to the x axis that crosses our current point and touches the smallest ellipse among the points the line. Similarly, when we update b we search the line that is parallel the bạz axis.
Following the points in Figure 4, one can see that the points obtained through our iterations in the one-variable case form similar triangles, where the ratio of corresponding sides for every two consecutive triangles are the same. Therefore, if we denote the side length of the first triangle with 1, then the side length for the next triangles are at az a. ,.... Using the g(2)
The sum over this geometric series is equal to I-a which inspires adaptive step size in notation in Figure 4, g0y even iterations.
Note that in the one-variable case, we can recover the optimal solution with only two iterations as the lines connecting (2,( in ((==))) go through the optimal solution. While this does not necessarily hold in higher-dimensional settings, experiments demonstrate that acceleration improves the number of iterations and the running time of our approach significantly, especially when the number of observed entries are small (i.e., when 8 is large).
18


Reviews:
[{"review_id": "IMonYi1sDD", "reviewer": "Reviewer_rkMi", "summary": "This paper introduces an efficient tensor completion algorithm, approximate-mini-ALS, which transforms unstructured tensor completion into a structured tensor decomposition problem using a lifting strategy and approximate Richardson iteration. By leveraging leverage score sampling, the method achieves sublinear time complexity and linear convergence with a contraction factor of \\(1 - 1/\\beta\\), theoretically guaranteeing near-optimal sampling complexity \\(O(n^{3/2})\\) ().", "questions": "", "limitations": "", "rating": 3, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your review and for recognizing that our approach is an innovative combination. Please see our responses to your questions and weaknesses below.\n\n> In the experimental section, the author presents too few experimental results, which do not provide sufficient evidence to convincingly demonstrate the advantages of the proposed algorithm. The author should refer to article and increase the number of experiments.\n\nWe believe our experiments provide sufficient evidence for the speed-up of our algorithms in practice (see the running time subplots in Figures 2 and 3). Please note that due to the improved asymptotic running time of our algorithm (several theorems in Section 4), the speed-up of our algorithm (relative to baselines) only gets faster for larger tensors. That said, we appreciate the suggestions for extra experiments and are curious which “article” you are referring to.\n\n> The author should provide ablation experiments\n\nWe believe our experiments include an ablation study, as we have provided experiments (1) with and without acceleration (i.e., the fourth column vs the third column in Figure 2; `accelerated-mini-als` is faster than `mini-als`), and (2) with and without leverage score sampling (i.e., `approximate-mini-als` vs `mini-als` in Figure 1). We appreciate more specific suggestions if there are ablations that you think would improve the paper.\n\n> Whether the method can be transferred to t-svd framework？\n\nSome recent studies have considered more efficient algorithms for t-SVD using leverage score sampling—see [1]. Due to this connection, we believe it might be possible to adopt our method for tensor completion under the t-SVD framework as well. We will add a more detailed discussion about this to the next version of the paper.\n\n\n[1] Tarzanagh, Davoud Ataee, and George Michailidis. \"Fast randomized algorithms for t-product based tensor operations and decompositions with applications to imaging data.\" SIAM Journal on Imaging Sciences 11, no. 4 (2018): 2629-2664.", "claims_and_evidence": "Please refer to below", "methods_and_evaluation_criteria": "Please refer to below", "theoretical_claims": "Please refer to below", "experimental_designs_or_analyses": "Please refer to below", "relation_to_broader_scientific_literature": "Please refer to below", "essential_references_not_discussed": "Please refer to below", "other_strengths_and_weaknesses": "Strengths：\n1.Innovative Combination: The lifting strategy combined with approximate Richardson iteration and leverage score sampling offers a novel approach to tensor completion, leveraging existing techniques in a non-trivial way.\n 2. Compared with ALS algorithms, the algorithm speed is improved significantly.\nWeaknesses:\n1.In the experimental section, the author presents too few experimental results, which do not provide sufficient evidence to convincingly demonstrate the advantages of the proposed algorithm. The author should refer to article and increase the number of experiments.", "other_comments_or_suggestions": "There are some typos，eg: mode- -72 unfolding...\nThe author should check typos carefully.", "questions_for_authors": "1.The author should provide ablation experiments.\n2.Whether the method can be transferred to t-svd framework？", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "bD4yXfg9Z4", "reviewer": "Reviewer_nT2U", "summary": "In this paper the well-known tensor completion problem for several popular low-rank tensor formats is considered. This problem is solved through a modified (randomized) ALS method, where it is assumed that we have access to all elements of the tensor (the so-called “lifting approach”). In those points whose values are not initially given, the values taken by the tensor at a given iteration is taken. Formally, the problem is mathematically written as a minimization problem not only on the elements of the expansion cores, but also on the unknown values of the tensor (the equivalence of these formulations is proved). One of the theoretical results of the paper is the construction of analogies of the obtained ALS iterations with Approximate Richardson Iteration. This allows authors to prove convergence theorems of the method. Small tensors (10^6~10^7 elements) are taken as numerical experiments and only the approximation time is compared with another methods (the accuracy of all methods is approximately the same).", "questions": "", "limitations": "", "rating": 3, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your review and for going over the theoretical results and experiments carefully. Please see our responses to your questions and weaknesses below.\n\n> The experiments presented are rather synthetic and small, it is not clear where in real-world applications these techniques can be applied\n\nNote that the `cardiac-mri` and `hyperspectral` tensors in Section 5.2 are real-world tensors, not synthetic.\n\n> Low-rank tensor decompositions were created to, among other things, overcome the “curse of dimensinality” --- situations where the full tensor does not fit in the memory of a computing system (or it is so large that it is difficult to work with). When the elements of the tensor occupy a small memory (in the experiments in the paper --- at most tens of megabytes) the need for low-rank decompositions as such is small; and the problem of restoring missing elements can be solved in many other ways (e.g., spline interpolation).\n\nWe completely agree that this is a primary motivation for low-rank tensor decompositions. That said, another benefit of low-rank decomposition is to prevent overfitting because it can also be viewed as a form of regularization.\n\nIn our experiments, our focus has been on comparison with baselines in terms of the running time and convergence. Our experiments demonstrate that we achieve a significant speed-up even for small tensors, and this speed-up only increases for larger tensors. Moreover, the memory footprint of our method is much smaller than the baselines due to the use of approximate solvers (using sampling techniques). In previous works on tensor decompositions that use leverage score sampling, it has been demonstrated that these methods are *especially effective* when direct methods run out of memory. Thus, we believe our experiments provide sufficient evidence for practicality and improvements of our algorithms.\n\n> For the same reason, the case when a very small fraction of elements of all possible tensor elements is known is interesting. In many real-world datasets this can be tenths and hundredths of a percent). And in this mode the presented algorithm performs poorly, as it can be seen from the plots.\n\nWe disagree with this assessment. Although the vanilla version of our algorithm is not faster than direct ALS methods when a small fraction of elements are observed, the *accelerated version* gives significant speed-ups even for this case—see the 4th column (`accelerated-mini-als`) of Figure 2 and compare it with the 3rd column (`mini-als`). Our asymptotic running times (theorems in Section 4) imply that such speed-ups will be even faster for larger tensors. In other words, for sufficiently large tensors with a small fraction of observed entries, our algorithm still provides significant speed-ups.\n\n> But for different algorithms, this time can strongly depend on the implementation, the use of hardware such as GPUs, the ability to parallelize and vectorize algorithms, etc. Thus, the only superiority of the presented algorithm in the experimental data, i.e., running time, is rather doubtful.\n\nOur algorithms similarly benefit from the use of parallelizability, vectorized algorithms, and the use of GPUs since they’re iterative methods under the hood. Note that `numpy.linalg` (and hence `tensorly`) calls BLAS and LAPACK low-level operations, some of which are parallelized on CPU. Further, improvements in asymptotic running times of algorithms often lead to speed-ups that might not be achievable with a better code or hardware. Our approach unlocks such improvements for tensor completion via techniques that have been developed for tensor decomposition.\n\n> Have you tried your method on much larger datasets (including those for which the full tensor does not fit into computer memory)?\n\nWe have not tried it on tensor datasets that exceed our machine’s memory. We believe this is a fantastic future direction of research, but it is somewhat out of the scope for this paper since such implementations require special care with respect to the code and hardware.\n\n> which parafac implementation did you use in experiments in Fig.2?\n\nPARAFAC refers to the EM approach of Tomasi-Bro (2005), which is equivalent to running our mini-ALS algorithm for exactly one step at each iteration of the ALS algorithm.", "claims_and_evidence": "The paper is largely theoretical, and it does have theorems for the convergence of the method and its running time. However, the experiments presented are rather synthetic and small, it is not clear where in real-world applications these techniques can be applied.", "methods_and_evaluation_criteria": "As numerical experiments, the authors take rather synthetic problems with tensors with a small number of elements. At the same time, the presented algorithm works well only when the ratio of the number of known elements to the number of unknown elements is rather large. \nDisadvantages of this approach:\n- Low-rank tensor decompositions were created to, among other things, overcome the “curse of dimensinality” --- situations where the full tensor does not fit in the memory of a computing system (or it is so large that it is difficult to work with).  When the elements of the tensor occupy a small memory (in the experiments in the paper --- at most tens of megabytes) the need for low-rank decompositions as such is small; and the problem of restoring missing elements can be solved in many other ways (e.g., spline interpolation).\n- For the same reason, the case when a very small fraction of elements of all possible tensor elements is known is interesting.  In many real-world datasets this can be tenths and hundredths of a percent). And in this mode the presented algorithm performs poorly, as it can be seen from the plots.\n- as a result, the authors show only the execution time of the algorithms (all of them have approximately the same accuracy). But for different algorithms, this time can strongly depend on the implementation, the use of hardware such as GPUs, the ability to parallelize and vectorize algorithms, etc. Thus, the only superiority of the presented algorithm in the experimental data, i.e., running time, is rather doubtful.", "theoretical_claims": "All theorems are proved correct in my opinion, theoretical results are true.", "experimental_designs_or_analyses": "Apart from the fundamental limitations of the experiments I wrote above, the numerical experiments in the paper are carried out without question.", "relation_to_broader_scientific_literature": "In fact, this paper is a some amalgamation of several existing and developed ideas --- 1. The ALS method, which can be applied to all three tensor decompositions presented (and in each of them it can be written as a regression problem), 2. Score Sampling, see e.g. the work of , Bharadwaj et al. (2024,) 3. Richardson Iterations. But ussually Score Sampling is applied in the case where any element of the tensor is accessible, not to the tensor competion problem.", "essential_references_not_discussed": "no", "other_strengths_and_weaknesses": "Minor:\n- l119: $:=R/R_k$ should not be in subscript", "other_comments_or_suggestions": "-", "questions_for_authors": "- Have you tried your method on much larger datasets (including those for which the full tensor does not fit into computer memory)?\n\n- which parafac implementation did you use in experiments in Fig.2?", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "z5BCQRonqV", "reviewer": "Reviewer_BWjF", "summary": "This paper addresses the problem of tensor completion (TC) by proposing a novel approach that leverages approximate Richardson iteration and structured tensor decomposition (TD) algorithms. The authors introduce a lifting technique to transform the TC problem into a structured linear regression problem, enabling the use of fast TD algorithms as blackbox subroutines. This approach allows for sublinear-time methods, significantly speeding up the completion process compared to traditional direct methods. The proposed method, termed approximate-mini-ALS, is theoretically analyzed for convergence and empirically validated on real-world tensors, demonstrating substantial speed improvements while maintaining comparable reconstruction accuracy.", "questions": "", "limitations": "", "rating": 4, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your review, commenting on its novelty, and carefully checking the proofs of the main lemmas/theorems. Please see our responses to your questions and weaknesses below.\n\n> I’m curious—when the number of observed values is small, does the author’s method have an advantage over traditional algorithms? Is it more efficient in terms of time, or does it have better recovery performance?\n\nPlease see the answer to the next question.\n\n> However, a weakness is that it remains unclear whether the proposed algorithm achieves faster speed or better recovery performance under high missing rate scenarios.\n\nThe vanilla version of our algorithm can be slower than direct ALS methods when a small fraction of elements are observed, but our accelerated version gives significant speed-ups even for this case—see the 4th column of Figure 2. Note that this occurs when $\\beta$ is large in Theorem 3.7, which happens if the lifted TD design matrix is not a tight spectral approximation to the TC design matrix. Our asymptotic running times imply that the speed-ups become even faster for larger tensors. In other words, for sufficiently large tensors with a small fraction of observed entries, our algorithm still provides significant speed-ups.\n\n> Additionally, the author mentions tensor networks in the supplementary material—could their method be applied to tensor network algorithms, such as FCTN decomposition or TW decomposition?\n\nYes, we believe so. Recently, [1] has proposed a framework for designing sampling-based decomposition algorithms for arbitrary tensor networks. Therefore, we think someone could use that approach to design more efficient algorithms for FCTN decomposition and then, using our methods, extend it to tensor completion. Even more recently, sampling-based algorithms were designed for TW decompositions [2]. Our approach can extend their algorithm to tensor completion as well. We will add a discussion about these in the next version.\n\n> x^{(k)} should represent the value from the k-th iteration of the ALS algorithm, but this is not clarified in the text. I believe this could confuse the readers, so please provide an explanation.\n\nThank you for the suggestion. We clarify this in the next version of the paper.\n\n> The author should provide the parameter settings used in the experiments or release the code to help readers better understand the work.\n\nNote that we included our code in the supplementary material. It is also available in a private GitHub repository. We will make the repository public and provide a link to it in the next version of the paper after the double-blind review process concludes.\n\n[1] Malik, Osman Asif, Vivek Bharadwaj, and Riley Murray. \"Sampling-based decomposition algorithms for arbitrary tensor networks.\" arXiv preprint arXiv:2210.03828 (2022).\n\n[2] Wang, Mengyu, Yajie Yu, and Hanyu Li. \"Randomized Tensor Wheel Decomposition.\" SIAM Journal on Scientific Computing 46, no. 3 (2024): A1714-A1746.", "claims_and_evidence": "The authors claim that their proposed method, approximate-mini-ALS, can achieve significant speedups in tensor completion tasks without compromising on the quality of the solution. The evidence supporting these claims includes:\n1.\tTheoretical Analysis: The authors provide a detailed convergence analysis of their approximate Richardson iteration-based algorithm, proving that it converges at the same rate as the exact Richardson iteration under certain conditions (Theorem 3.7).\n2.\tEmpirical Validation: Extensive experiments on synthetic and real-world tensors (e.g., CARDIAC-MRI and HYPERSPECTRAL datasets) demonstrate that the proposed method can be orders of magnitude faster than direct methods while achieving similar reconstruction errors.\n3.\tComparison with Existing Methods: The paper compares the proposed method with existing approaches such as direct ALS and the EM-based algorithm of Tomasi & Bro (2005), showing superior performance in terms of running time and solution quality.", "methods_and_evaluation_criteria": "The proposed method-approximate-mini-ALS using RRE to measure the relative Frobenius norm of the difference between the completed tensor and the ground truth.I think there is no problem.", "theoretical_claims": "Yes, I check the proof of Lemma 3.1, 3.2, 3.5, Theorem 3.7. They are all correct.", "experimental_designs_or_analyses": "I think experiments are reasonably well designed. In addition, TC problem is often used for predicting missing values. The author states that the complexity of other algorithms is related to the number of observed values. I’m curious—when the number of observed values is small, does the author’s method have an advantage over traditional algorithms? Is it more efficient in terms of time, or does it have better recovery performance?", "relation_to_broader_scientific_literature": "The paper builds upon and contributes to several areas of the broader scientific literature:\n1.\tTensor Completion: Extends the state-of-the-art in TC by proposing a novel lifting approach that leverages structured TD algorithms, addressing the limitations of existing methods that rely on direct matrix operations.\n2.\tIterative Methods: Connects TC to the Richardson iteration, a classical iterative method for solving linear systems, providing new insights into the convergence properties of TC algorithms.\n3.\tTensor Decompositions: Integrates with existing TD methods (e.g., CP, Tucker, tensor-train) by using them as subroutines, ensuring compatibility with a wide range of tensor structures and applications.", "essential_references_not_discussed": "No", "other_strengths_and_weaknesses": "The strength of this paper is that it proposes a novel and fast tensor completion algorithm with originality, and provides theoretical proofs to support it. Experimental results demonstrate its effectiveness. However, a weakness is that it remains unclear whether the proposed algorithm achieves faster speed or better recovery performance under high missing rate scenarios. Additionally, the author mentions tensor networks in the supplementary material—could their method be applied to tensor network algorithms, such as FCTN decomposition or TW decomposition?", "other_comments_or_suggestions": "•  x^{(k)} should represent the value from the k-th iteration of the ALS algorithm, but this is not clarified in the text. I believe this could confuse the readers, so please provide an explanation.\n•  The author should provide the parameter settings used in the experiments or release the code to help readers better understand the work.", "questions_for_authors": "1. The author should provide the parameter settings used in the experiments or release the code to help readers better understand the work.\n2. TC problem is often used for predicting missing values. The author states that the complexity of other algorithms is related to the number of observed values. I’m curious—when the number of observed values is small, does the author’s method have an advantage over traditional algorithms? Is it more efficient in terms of time, or does it have better recovery performance?\n3. The author mentions tensor networks in the supplementary material—could their method be applied to tensor network algorithms, such as FCTN decomposition or TW decomposition?", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "2utBWaTAmr", "reviewer": "Reviewer_NbQN", "summary": "This paper addresses the computational challenges of Tensor Completion (TC) by proposing a novel method that integrates Tensor Decompositions (TDs)—including CANDECOMP/PARAFAC (CP), Tucker, and Tensor Train (TT)—with a lifting framework. \n\nThe authors reformulate TC as a structured tensor decomposition problem, effectively reducing its computational complexity. This reformulation enables the application of efficient Alternating Least Squares (ALS) algorithms. By interpreting the algorithm as a preconditioned Richardson iteration, the authors establish rigorous theoretical convergence guarantees. \n\nTo further enhance efficiency, the method incorporates sketching techniques via leverage score sampling, strategically subsampling tensor fibers to accelerate computations. \n\nComprehensive experiments on synthetic and real-world datasets demonstrate that the proposed approach achieves significant speedups over state-of-the-art TC methods while maintaining competitive accuracy. The integration of lifting with sketching and convergence analysis presents a key advancement in scalable tensor completion.\n\n## update after rebuttal\nI thank very much the authors for clarifying the issues I pointed out and their promise to improve the final version accordingly. I keep my original score.", "questions": "", "limitations": "", "rating": 4, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for your review and careful read. Please see our responses to your questions and weaknesses below.\n\n> The general TC problem of eq. (1) is not alligned with the approach of the current paper where a TD approach is used meaning that the approximation of the data tensor is a low-rank tensor decomposition (CP, Tucker or TT) where its factors are the minimizing parameter \\theta. In other words, the regularization term of eq. (1) is not needed nor used in this paper.\n\nWe agree that in the algorithm explicitly discussed in the paper, the regularization term of eq. (1) is not used. However, we emphasize that our approach can be applied to methods with regularization as well if the corresponding regression problems in such methods are structured and that structure enables faster algorithms. One example is Tucker decomposition with $\\ell_2$-regularization, which reduces to solving *ridge regression* problems. For this problem, fast algorithms using sampling techniques are developed in [1], and our approach gives speed-ups for the corresponding tensor completion problem. Thank you for pointing this out — we will add a discussion about this to the next version.\n\n> The statement “Rank constraints can be incorporated into (1) by including appropriate penalty terms in R(\\theta)” is confusing since in this paper, the low rank constraint is explicitly imposed by chosing the rank of the decomposition.\n\nWe completely agree that the algorithms we explicitly describe in the paper do not incorporate rank constraints this way. The comment is meant to state that if there are regularization terms that still lead to structured problems that can be solved faster, then our lifting approach can give speed-ups for them as well, e.g., ridge regression that is discussed in the answer to the previous question. Note that the ridge regression reduces the effective dimension of the problem and can be used to implicitly induce a lower rank for the decomposition. We will make this clearer in the next version to avoid confusion.\n\n> In the paper, it is not clearly explained what the authors refers as direct and parafac methods.\n\nThis is explained in the first column on page 8. The `direct` method uses normal equations for solving the linear regression problem in each iteration of the ALS algorithm. More specifically for a matrix $A$, vector $b$, and set of observed entries $\\Omega$, we require to solve $x^*=\\arg\\min_{x} || A_{\\Omega} x - b_{\\Omega} ||_2^2$. \n\nThe `direct` method computes $x^*$ using the following operation $x^* = (A_{\\Omega}^{\\top} A_{\\Omega})^{-1} A_{\\Omega}^{\\top} b_{\\Omega}$. \n\nThe `parafac` method is the EM approach of Tomasi-Bro (2005), which is equivalent to running the mini-ALS algorithm for one step at each iteration of the ALS algorithm. We will add more detailed explanations about these algorithms in the next version of the paper.\n\n> What is the difference between the plots in 3rd and 4th columns in Figure 2? They both compare the running time of the evaluated algorithms as a function of sample ratio. Do they refer to mini-ALS and accelerated-mini-ALS, respectively? Please, confirm.\n\nYes, the third column in Figure 2 shows the running times of the `mini-ALS` algorithm for different values of $\\varepsilon$, while the fourth column shows the running times of the `approximate-mini-ALS` algorithm. This is explained in the text (lines 407–425, second column on page 8). We will add a comprehensive explanation to the caption of Figure 2 in the next version.\n\n> Figure 1 show that the approximate-mini-als algorithm achieves lower error than the other algorithms for a range of steps around step 15. This is surprising and should be explained.\n\nThe `approximate-mini-ALS` in Figure 1 is a randomized method that uses leverage score sampling for the Kronecker product. We believe the lower error is due to randomness, which causes the algorithm to follow a different convergence path. We will provide an explanation of this in the next version of the paper.\n\nReferences:\n\n[1] Fahrbach, Matthew, Gang Fu, and Mehrdad Ghadiri. \"Subquadratic kronecker regression with applications to tensor decomposition.\" Advances in Neural Information Processing Systems 35 (2022): 28776-28789.", "claims_and_evidence": "The claims in the submission are generally well-supported by theoretical analysis and empirical results. The main claims are:\n\n-\tReduction in Computational Complexity: The paper provides theoretical justifications for the efficiency of the proposed method using leverage score sampling and the lifting approach. Running time analyses for CP, Tucker, and TT decompositions are given. Empirical results show significant speed improvements over standard ALS-based CP TC methods.\n\n-\tInterpretation as a Preconditioned Richardson Iteration: Theoretical derivations (Lemma 3.5) show that the mini-ALS algorithm simulates Richardson iteration.", "methods_and_evaluation_criteria": "The proposed methods and evaluation criteria are reasonable and relevant for the tensor completion (TC) problem. However, certain aspects require further clarification (see detailed comments below).", "theoretical_claims": "I have not verified the proofs provided in the Appendices. However, the theoretical claims appear well-founded.", "experimental_designs_or_analyses": "The experimental designs and analysis are correct and useful. However, further clarification is needed (see detailed comments below).", "relation_to_broader_scientific_literature": "The paper effectively integrates concepts from tensor decomposition, iterative optimization, and randomized numerical linear algebra. By leveraging preconditioned Richardson iteration and leverage score sampling, it advances the field of efficient tensor completion while maintaining theoretical guarantees. \n\nHowever, it focuses solely on tensor decomposition-based completion, while many other tensor completion methods exist, including deep learning-based approaches and nuclear norm minimization, among others.", "essential_references_not_discussed": "No, I think the papers included in the references are a good selection of relevant works in this field.", "other_strengths_and_weaknesses": "Strengths:\n\n-\tTheoretical Contribution: A useful connection between TC and preconditioned Richardson iteration, providing rigorous convergence guarantees (Theorem 3.7) is provided. The lifting approach restores structure in TC regression problems, enabling the use of fast TD algorithms as subroutines. This bridges iterative optimization and tensor algebra in a principled manner.\n\n-\tGenerality: The framework applies to multiple tensor decompositions (CP, Tucker, TT), demonstrating flexibility. By using TD algorithms as black-box subroutines, the method benefits from future advancements in TD research.\n\n-\tPractical Efficiency: Experiments on real-world datasets (e.g., cardiac MRI, hyperspectral imaging) show significant speedups (up to 100× faster than direct methods) while maintaining competitive reconstruction errors.\n\nWeaknesses:\n\n-\tClarity and Accessibility: some parts requires further clarification (see my detailed comments in the next section)\n\n-\tComparison to State-of-the-Art: the paper does not benchmark against modern TC methods (e.g., tensor nuclear norm minimization, deep learning-based approaches), leaving its relative performance unclear.", "other_comments_or_suggestions": "Here, I provide my detailed comments/suggestions to improve clarity and correctness of the paper:\n\n-\tThe general TC problem of eq. (1) is not alligned with the approach of the current paper where a TD approach is used meaning that the approximation of the data tensor is a low-rank tensor decomposition (CP, Tucker or TT) where its factors are the minimizing parameter \\theta. In other words, the regularization term of eq. (1) is not needed nor used in this paper. \n\n-\tThe statement “Rank constraints can be incorporated into (1) by including appropriate penalty terms in R(\\theta)” is confusing since in this paper, the low rank constraint is explicitly imposed by chosing the rank of the decomposition.\n\n-\tExperimental results requires further clarification:\n\n1)\tFigure 1 show that the approximate-mini-als algorithm achieves lower error than the other algorithms for a range of steps around step 15. This is surprising and should be explained.\n\n2)\tIn the paper, it is not clearly explained what the authors refers as direct and parafac methods.\n\n3)\tFigure 2 is confusing, for example there is no explanation about what do the plots in 3rd and 4th columns refer to. The caption should be improved.", "questions_for_authors": "I have some questions that could help me to better understand the experimental results:\n\n-\tWhy  does approximate-mini-als achieve lower error than other algorithms around step 15 in Fig. 1?\n\n-\tCould you please explain what do you mean the direct and parafac methods used in Fig. 2? You should precisely define those methods.\n\n-\tWhat is the difference between the plots in 3rd and 4th columns in Figure 2? They both compare the running time of the evaluated algorithms as a function of sample ratio. Do they refer to mini-ALS and accelerated-mini-ALS, respectively? Please, confirm.", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Meta-review: The paper reviewers all agree that the paper is worth publishing at ICML. The combination of methods in the paper seems innovative and the paper seems like a good contribution.
Decision: Accept (poster)