Thus, applying Theorem C.1 by fixing $U^{\prime},V^{\prime},W^{\prime}$ and
using Theorem C.45 at the end, we can output three matrices $U,V,W$, where
each entry can be written using $n^{\delta^{\prime}}$ bits, so that we satisfy
$\|U\otimes V\otimes W-A\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}$.

For the running time, since each entry of the input is bounded by $n^{\delta}$
bits, due to Theorem C.1, we need
$(\operatorname{nnz}(A)+n\operatorname{poly}(k/\epsilon))\cdot n^{\delta}$
time to reduce the size of the problem to $\operatorname{poly}(k/\epsilon)$
size (with each number represented using $O(n^{\delta})$ bits). According to
Theorem C.45, the running time of using a polynomial system verifier to get
the solution is
$2^{{O}(k^{2}/\epsilon)}n^{O(\delta^{\prime})}=2^{{O}(k^{2}/\epsilon)}n^{O(\delta)}$
time. Thus the total running time is
$(\operatorname{nnz}(A)+n\operatorname{poly}(k/\epsilon))n^{\delta}+2^{{O}(k^{2}/\epsilon)}\cdot
n^{O(\delta)}$.

Part (ii@) is similar to Part (i@). Part (iii@) is trivial to prove since
there exists a solution which can be written down in the bit model, so we
obtain a $(1+\epsilon)$-approximation. Part (iv@) is also very similar to Part
(ii@).

∎

###  C.2 Input sparsity reduction

Algorithm 3 Reducing the Size of the Objective Function from
$\operatorname{poly}(n)$ to $\operatorname{poly}(k)$

1:procedure
FInputSparsityReduction($A,V_{1},V_{2},V_{3},n,b_{1},b_{2},b_{3},k,\epsilon$)
$\triangleright$ Lemma C.3

2: $c_{1}\leftarrow c_{2}\leftarrow
c_{3}\leftarrow\operatorname{poly}(k,1/\epsilon)$.

3: Choose sparse embedding matrices $T_{1}\in\mathbb{R}^{c_{1}\times n}$,
$T_{2}\in\mathbb{R}^{c_{2}\times n}$, $T_{3}\in\mathbb{R}^{c_{3}\times n}$.
$\triangleright$ Definition B.16

4: $\widehat{V}_{i}\leftarrow T_{i}V_{i}\in\mathbb{R}^{c_{i}\times
b_{i}},\forall i\in[3]$.

5: $C\leftarrow A(T_{1},T_{2},T_{3})\in\mathbb{R}^{c_{1}\times c_{2}\times
c_{3}}$.

6: return $\widehat{V}_{1}$, $\widehat{V}_{2}$, $\widehat{V}_{3}$ and $C$.

7:end procedure

######  Lemma C.3.

Let $\operatorname{poly}(k,1/\epsilon)\geq b_{1}b_{2}b_{3}\geq k$. Given a
tensor $A\in\mathbb{R}^{n\times n\times n}$ and three matrices
$V_{1}\in\mathbb{R}^{n\times b_{1}}$, $V_{2}\in\mathbb{R}^{n\times b_{2}}$,
and $V_{3}\in\mathbb{R}^{n\times b_{3}}$, there exists an algorithm that takes
$O(\operatorname{nnz}(A)+\operatorname{nnz}(V_{1})+\operatorname{nnz}(V_{2})+\operatorname{nnz}(V_{3}))=O(\operatorname{nnz}(A)+n\operatorname{poly}(k/\epsilon))$
time and outputs a tensor $C\in\mathbb{R}^{c_{1}\times c_{2}\times c_{3}}$ and
three matrices $\widehat{V}_{1}\in\mathbb{R}^{c_{1}\times b_{1}}$,
$\widehat{V}_{2}\in\mathbb{R}^{c_{2}\times b_{2}}$ and
$\widehat{V}_{3}\in\mathbb{R}^{c_{3}\times b_{3}}$ with
$c_{1}=c_{2}=c_{3}=\operatorname{poly}(k,1/\epsilon)$, such that with
probability at least $0.99$, for all
$\alpha>0,X_{1},X^{\prime}_{1}\in\mathbb{R}^{b_{1}\times
k},X_{2},X^{\prime}_{2}\in\mathbb{R}^{b_{2}\times
k},X_{3},X^{\prime}_{3}\in\mathbb{R}^{b_{3}\times k}$ satisfy that,

$\displaystyle\left\|\sum_{i=1}^{k}(\widehat{V}_{1}X_{1}^{\prime})_{i}\otimes(\widehat{V}_{2}X_{2}^{\prime})_{i}\otimes(\widehat{V}_{3}X_{3}^{\prime})_{i}-C\right\|_{F}^{2}\leq\alpha\left\|\sum_{i=1}^{k}(\widehat{V}_{1}X_{1})_{i}\otimes(\widehat{V}_{2}X_{2})_{i}\otimes(\widehat{V}_{3}X_{3})_{i}-C\right\|_{F}^{2},$

then,

$\displaystyle\left\|\sum_{i=1}^{k}(V_{1}X_{1}^{\prime})_{i}\otimes(V_{2}X_{2}^{\prime})_{i}\otimes(V_{3}X_{3}^{\prime})_{i}-A\right\|_{F}^{2}\leq(1+\epsilon)\alpha\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

###### Proof.

Let $X_{1}\in\mathbb{R}^{b_{1}\times k},X_{2}\in\mathbb{R}^{b_{2}\times
k},X_{3}\in\mathbb{R}^{b_{3}\times k}.$ First, we define
$Z_{1}=((V_{2}X_{2})^{\top}\odot(V_{3}X_{3})^{\top})\in\mathbb{R}^{k\times
n^{2}}$. (Note that, for each $i\in[k]$, the $i$-th row of matrix $Z_{1}$ is
$\operatorname{vec}((V_{2}X_{2})_{i}\otimes(V_{3}X_{3})_{i})$.) Then, by
flattening we have

$\displaystyle\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}=\|V_{1}X_{1}\cdot
Z_{1}-A_{1}\|_{F}^{2}.$

We choose a sparse embedding matrix (Definition B.16)
$T_{1}\in\mathbb{R}^{c_{1}\times n}$ with
$c_{1}=\operatorname{poly}(k,1/\epsilon)$ rows. Since $V_{1}$ has
$b_{1}\leq\operatorname{poly}(k/\epsilon)$ columns, according to Lemma B.19
with probability $0.999$, for all $X_{1}\in\mathbb{R}^{b_{1}\times
k},Z\in\mathbb{R}^{k\times n^{2}}$,

$\displaystyle(1-\epsilon)\|V_{1}X_{1}Z-A_{1}\|_{F}^{2}\leq\|T_{1}V_{1}X_{1}Z-T_{1}A_{1}\|_{F}^{2}\leq(1+\epsilon)\|V_{1}X_{1}Z-A_{1}\|_{F}^{2}.$

Therefore, we have

$\displaystyle\|T_{1}V_{1}X_{1}\cdot
Z_{1}-T_{1}A_{1}\|_{F}^{2}=(1\pm\epsilon)\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

Second, we unflatten matrix $T_{1}A_{1}\in\mathbb{R}^{c_{1}\times n^{2}}$ to
obtain a tensor $A^{\prime}\in\mathbb{R}^{c_{1}\times n\times n}$. Then we
flatten $A^{\prime}$ along the second direction to obtain
$A_{2}\in\mathbb{R}^{n\times c_{1}n}$. We define
$Z_{2}=(T_{1}V_{1}X_{1})^{\top}\odot(V_{3}X_{3})^{\top}\in\mathbb{R}^{k\times
c_{1}n}$. Then, by flattening,

$\displaystyle\|V_{2}X_{2}\cdot Z_{2}-A_{2}\|_{F}^{2}=$
$\displaystyle\leavevmode\nobreak\ \|T_{1}V_{1}X_{1}\cdot
Z_{1}-T_{1}A_{1}\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
(1\pm\epsilon)\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

We choose a sparse embedding matrix (Definition B.16)
$T_{2}\in\mathbb{R}^{c_{2}\times n}$ with
$c_{2}=\operatorname{poly}(k,1/\epsilon)$ rows. Then according to Lemma B.19
with probability $0.999$, for all $X_{2}\in\mathbb{R}^{b_{2}\times k}$,
$Z\in\mathbb{R}^{k\times c_{1}n}$,

$\displaystyle(1-\epsilon)\|V_{2}X_{2}Z-A_{2}\|_{F}^{2}\leq\|T_{2}V_{2}X_{2}Z-T_{2}A_{2}\|_{F}^{2}\leq(1+\epsilon)\|V_{2}X_{2}Z-A_{2}\|_{F}^{2}.$

Therefore, we have

$\displaystyle\|T_{2}V_{2}X_{2}\cdot Z_{2}-T_{2}A_{2}\|_{F}^{2}=$
$\displaystyle\leavevmode\nobreak\ (1\pm\epsilon)\|V_{2}X_{2}\cdot
Z_{2}-A_{2}\|_{F}^{2}$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
(1\pm\epsilon)^{2}\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

Third, we unflatten matrix $T_{2}A_{2}\in\mathbb{R}^{c_{2}\times c_{1}n}$ to
obtain a tensor $A^{\prime\prime}(=A(T_{1},T_{2},I))\in\mathbb{R}^{c_{1}\times
c_{2}\times n}$. Then we flatten tensor $A^{\prime\prime}$ along the last
direction (the third direction) to obtain matrix $A_{3}\in\mathbb{R}^{n\times
c_{1}c_{2}}$. We define
$Z_{3}=(T_{1}V_{1}X_{1})^{\top}\odot(T_{2}V_{2}X_{2})^{\top}\in\mathbb{R}^{k\times
c_{1}c_{2}}$. Then, by flattening, we have

$\displaystyle\|V_{3}X_{3}\cdot Z_{3}-A_{3}\|_{F}^{2}=$
$\displaystyle\leavevmode\nobreak\ \|T_{2}V_{2}X_{2}\cdot
Z_{2}-T_{2}A_{2}\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
(1\pm\epsilon)^{2}\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

We choose a sparse embedding matrix (Definition B.16)
$T_{3}\in\mathbb{R}^{c_{3}\times n}$ with
$c_{3}=\operatorname{poly}(k,1/\epsilon)$ rows. Then according to Lemma B.19
with probability $0.999$, for all $X_{3}\in\mathbb{R}^{b_{3}\times k}$,
$Z\in\mathbb{R}^{k\times c_{1}c_{2}}$,

$\displaystyle(1-\epsilon)\|V_{3}X_{3}Z-A_{3}\|_{F}^{2}\leq\|T_{3}V_{3}X_{3}Z-T_{3}A_{3}\|_{F}^{2}\leq(1+\epsilon)\|V_{3}X_{3}Z-A_{3}\|_{F}^{2}.$

Therefore, we have

$\displaystyle\|T_{3}V_{3}X_{3}\cdot
Z_{3}-T_{3}A_{3}\|_{F}^{2}=(1\pm\epsilon)^{3}\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

Note that

$\displaystyle\|T_{3}V_{3}X_{3}\cdot
Z_{3}-T_{3}A_{3}\|_{F}^{2}=\left\|\sum_{i=1}^{k}(T_{1}V_{1}X_{1})_{i}\otimes(T_{2}V_{2}X_{2})_{i}\otimes(T_{3}V_{3}X_{3})_{i}-A(T_{1},T_{2},T_{3})\right\|_{F}^{2},$

and thus, we have $\forall X_{1}\in\mathbb{R}^{b_{1}\times
k},X_{2}\in\mathbb{R}^{b_{2}\times k},X_{3}\in\mathbb{R}^{b_{3}\times k}$

$\displaystyle\left\|\sum_{i=1}^{k}(T_{1}V_{1}X_{1})_{i}\otimes(T_{2}V_{2}X_{2})_{i}\otimes(T_{3}V_{3}X_{3})_{i}-A(T_{1},T_{2},T_{3})\right\|_{F}^{2}$
$\displaystyle=$
$\displaystyle(1\pm\epsilon)^{3}\left\|\sum_{i=1}^{k}({V}_{1}X_{1})_{i}\otimes({V}_{2}X_{2})_{i}\otimes({V}_{3}X_{3})_{i}-A\right\|_{F}^{2}.$

Let $\widehat{V}_{i}$ denote $T_{i}V_{i}$, for each $i\in[3]$. Let
$C\in\mathbb{R}^{c_{1}\times c_{2}\times c_{3}}$ denote
$A(T_{1},T_{2},T_{3})$. For $\alpha>1$, if

$\displaystyle\left\|\sum_{i=1}^{k}(\widehat{V}_{1}X_{1}^{\prime})_{i}\otimes(\widehat{V}_{2}X_{2}^{\prime})_{i}\otimes(\widehat{V}_{3}X_{3}^{\prime})_{i}-C\right\|_{F}^{2}\leq\alpha\left\|\sum_{i=1}^{k}(\widehat{V}_{1}X_{1})_{i}\otimes(\widehat{V}_{2}X_{2})_{i}\otimes(\widehat{V}_{3}X_{3})_{i}-C\right\|_{F}^{2},$

then

$\displaystyle\leavevmode\nobreak\
\left\|\sum_{i=1}^{k}(V_{1}X_{1}^{\prime})_{i}\otimes(V_{2}X_{2}^{\prime})_{i}\otimes(V_{3}X_{3}^{\prime})_{i}-C\right\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\frac{1}{(1-\epsilon)^{3}}\left\|\sum_{i=1}^{k}(\widehat{V}_{1}X_{1}^{\prime})_{i}\otimes(\widehat{V}_{2}X_{2}^{\prime})_{i}\otimes(\widehat{V}_{3}X_{3}^{\prime})_{i}-C\right\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\frac{1}{(1-\epsilon)^{3}}\alpha\left\|\sum_{i=1}^{k}(\widehat{V}_{1}X_{1})_{i}\otimes(\widehat{V}_{2}X_{2})_{i}\otimes(\widehat{V}_{3}X_{3})_{i}-C\right\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\frac{(1+\epsilon)^{3}}{(1-\epsilon)^{3}}\alpha\left\|\sum_{i=1}^{k}(V_{1}X_{1})_{i}\otimes(V_{2}X_{2})_{i}\otimes(V_{3}X_{3})_{i}-C\right\|_{F}^{2}$

By rescaling $\epsilon$ by a constant, we complete the proof of correctness.

##### Running time.

According to Section B.6, for each $i\in[3]$, $T_{i}V_{i}$ can be computed in
$O(\operatorname{nnz}(V_{i}))$ time, and $A(T_{1},T_{2},T_{3})$ can be
computed in $O(\operatorname{nnz}(A))$ time.

By the analysis above, the proof is complete. ∎

###  C.3 Tensor multiple regression

Algorithm 4 Frobenius Norm Tensor Multiple Regression

1:procedure FTensorMultipleRegression($A,U,V,d,n,k$) $\triangleright$ Theorem
C.4

2: $s\leftarrow O(k^{2}+k/\epsilon)$.

3: Choose $S\in\mathbb{R}^{n^{2}\times s}$ to be a TensorSketch.
$\triangleright$ Definition B.34

4: Compute $A\cdot S$.

5: Compute $B\cdot S$. $\triangleright$ $B=U^{\top}\odot V^{\top}$

6: $W\leftarrow(AS)(BS)^{\dagger}$

7: return $W$.

8:end procedure

######  Theorem C.4.

Given matrices $A\in\mathbb{R}^{d\times n^{2}}$, $U,V\in\mathbb{R}^{n\times
k}$, let $B\in\mathbb{R}^{k\times n^{2}}$ denote $U^{\top}\odot V^{\top}$.
There exists an algorithm that takes
$O(\operatorname{nnz}(A)+\operatorname{nnz}(U)+\operatorname{nnz}(V)+d\operatorname{poly}(k,1/\epsilon))$
time and outputs a matrix $W^{\prime}\in\mathbb{R}^{d\times k}$ such that,

$\displaystyle\|W^{\prime}B-A\|_{F}^{2}\leq(1+\epsilon)\min_{W\in\mathbb{R}^{d\times
k}}\|WB-A\|_{F}^{2}.$

###### Proof.

We choose a TensorSketch (Definition B.34) $S\in\mathbb{R}^{n^{2}\times s}$ to
reduce the problem to a smaller problem,

$\displaystyle\min_{W\in\mathbb{R}^{d\times k}}\|WBS-AS\|_{F}^{2}.$

Let $W^{\prime}$ denote the optimal solution to the above problem. Following a
similar proof to that in Section C.7.3, if $S$ is a $(1\pm 1/2)$-subspace
embedding and satisfies $\sqrt{\epsilon/k}$-approximate matrix product, then
$W^{\prime}$ provides a $(1+\epsilon)$-approximation to the original problem.
By Theorem B.35, we have $s=O(k^{2}+k/\epsilon)$.

##### Running time.

According to Definition B.34, $BS$ can be computed in
$O(\operatorname{nnz}(U)+\operatorname{nnz}(V))+\operatorname{poly}(k/\epsilon)$
time. Notice that each row of $S$ has exactly $1$ nonzero entry, thus $AS$ can
be computed in $O(\operatorname{nnz}(A))$ time. Since
$BS\in\mathbb{R}^{k\times s}$ and $AS\in\mathbb{R}^{d\times s}$,
$\min_{W\in\mathbb{R}^{d\times k}}\|WBS-AS\|_{F}^{2}$ can be solved in
$d\operatorname{poly}(sk)=d\operatorname{poly}(k/\epsilon)$ time. ∎

###  C.4 Bicriteria algorithms

####  C.4.1 Solving a small regression problem

######  Lemma C.5.

Given tensor $A\in\mathbb{R}^{n\times n\times n}$ and three matrices
$U\in\mathbb{R}^{n\times s_{1}},V\in\mathbb{R}^{n\times s_{2}}$ and
$W\in\mathbb{R}^{n\times s_{3}}$, there exists an algorithm that takes
$O(\operatorname{nnz}(A)+n\operatorname{poly}(s_{1},s_{2},s_{3},1/\epsilon))$
time and outputs $\alpha^{\prime}\in\mathbb{R}^{s_{1}\times s_{2}\times
s_{3}}$ such that

$\displaystyle\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha^{\prime}_{i,j,l}\cdot
U_{i}\otimes V_{j}\otimes
W_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\min_{\alpha\in\mathbb{R}^{s_{1}\times
s_{2}\times
s_{3}}}\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha_{i,j,l}\cdot
U_{i}\otimes V_{j}\otimes W_{l}-A\right\|_{F}^{2}.$

holds with probability at least $.99$.

###### Proof.

We define $\widetilde{b}\in\mathbb{R}^{n^{3}}$ to be the vector where the
$i+(j-1)n+(l-1)n^{2}$-th entry of $\widetilde{b}$ is $A_{i,j,l}$. We define
$\widetilde{A}\in\mathbb{R}^{n^{3}\times s_{1}s_{2}s_{3}}$ to be the matrix
where the
$(i+(j-1)n+(l-1)n^{2},i^{\prime}+(j^{\prime}-1)s_{2}+(l^{\prime}-1)s_{2}s_{3})$
entry is $U_{i^{\prime},i}\cdot V_{j^{\prime},j}\cdot W_{l^{\prime},l}$. This
problem is equivalent to a linear regression problem,

$\displaystyle\min_{x\in\mathbb{R}^{s_{1}s_{2}s_{3}}}\|\widetilde{A}x-\widetilde{b}\|_{2}^{2},$

where $\widetilde{A}\in\mathbb{R}^{n^{3}\times
s_{1}s_{2}s_{3}},\widetilde{b}\in\mathbb{R}^{n^{3}}$. Thus, it can be solved
fairly quickly using recent work [CW13, MM13, NN13]. However, the running time
of this naïvely is $\Omega(n^{3})$, since we have to write down each entry of
$\widetilde{A}$. In the next few paragraphs, we show how to improve the
running time to
$\operatorname{nnz}(A)+n\operatorname{poly}(s_{1},s_{2},s_{3})$.

Since $\alpha\in\mathbb{R}^{s_{1}\times s_{2}\times s_{3}},$ $\alpha$ can be
always written as $\alpha=X_{1}\otimes X_{2}\otimes X_{3},$ where
$X_{1}\in\mathbb{R}^{s_{1}\times
s_{1}s_{2}s_{3}},X_{2}\in\mathbb{R}^{s_{2}\times
s_{1}s_{2}s_{3}},X_{3}\in\mathbb{R}^{s_{3}\times s_{1}s_{2}s_{3}}$, we have

$\displaystyle\min_{\alpha\in\mathbb{R}^{s_{1}\times s_{2}\times
s_{3}}}\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha_{i,j,l}\cdot
U_{i}\otimes V_{j}\otimes
W_{l}-A\right\|_{F}^{2}=\underset{X_{3}\in\mathbb{R}^{s_{3}\times
s_{1}s_{2}s_{3}}}{\underset{X_{2}\in\mathbb{R}^{s_{2}\times
s_{1}s_{2}s_{3}}}{\underset{X_{1}\in\mathbb{R}^{s_{1}\times
s_{1}s_{2}s_{3}}}{\min}}}\left\|(UX_{1})\otimes(VX_{2})\otimes(WX_{3})-A\right\|_{F}^{2}.$

By Lemma C.3, we can reduce the problem size $n\times n\times n$ to a smaller
problem that has size $t_{1}\times t_{2}\times t_{3}$,

$\displaystyle\leavevmode\nobreak\
\min_{X_{1},X_{2},X_{3}}\left\|\sum_{i=1}^{s_{1}s_{2}s_{3}}(T_{1}UX_{1})_{i}\otimes(T_{2}VX_{2})_{i}\otimes(T_{3}WX_{3})_{i}-A(T_{1},T_{2},T_{3})\right\|_{F}^{2}$

where $T_{1}\in\mathbb{R}^{t_{1}\times n}$, $T_{2}\in\mathbb{R}^{t_{2}\times
n}$, $T_{3}\in\mathbb{R}^{t_{3}\times
n},t_{1}=t_{2}=t_{3}=\operatorname{poly}(s_{1}s_{2}s_{3}/\epsilon)$. Notice
that

$\displaystyle\leavevmode\nobreak\
\min_{X_{1},X_{2},X_{3}}\left\|\sum_{i=1}^{s_{1}s_{2}s_{3}}(T_{1}UX_{1})_{i}\otimes(T_{2}VX_{2})_{i}\otimes(T_{3}WX_{3})_{i}-A(T_{1},T_{2},T_{3})\right\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{\alpha\in\mathbb{R}^{s_{1}\times s_{2}\times
s_{3}}}\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha_{i,j,l}\cdot(T_{1}U)_{i}\otimes(T_{2}V)_{j}\otimes(T_{3}W)_{l}-A(T_{1},T_{2},T_{3})\right\|_{F}^{2}.$

Let

$\displaystyle\alpha^{\prime}=\underset{\alpha\in\mathbb{R}^{s_{1}\times
s_{2}\times
s_{3}}}{\arg\min}\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha_{i,j,l}\cdot(T_{1}U)_{i}\otimes(T_{2}V)_{j}\otimes(T_{3}W)_{l}-A(T_{1},T_{2},T_{3})\right\|_{F}^{2},$

then we have

$\displaystyle\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha^{\prime}_{i,j,l}\cdot
U_{i}\otimes V_{j}\otimes
W_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\min_{\alpha\in\mathbb{R}^{s_{1}\times
s_{2}\times
s_{3}}}\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}\sum_{l=1}^{s_{3}}\alpha_{i,j,l}\cdot
U_{i}\otimes V_{j}\otimes W_{l}-A\right\|_{F}^{2}.$

Again, according to Lemma C.3, the total running time is then
$O(\operatorname{nnz}(A)+n\operatorname{poly}(s_{1},s_{2},s_{3},1/\epsilon))$.

∎

######  Lemma C.6.

Given tensor $A\in\mathbb{R}^{n\times n\times n}$, and two matrices
$U\in\mathbb{R}^{n\times s},V\in\mathbb{R}^{n\times s}$ with
$\operatorname{rank}(U)=r_{1},\operatorname{rank}(V)=r_{2}$, let
$T_{1}\in\mathbb{R}^{t_{1}\times n},T_{2}\in\mathbb{R}^{t_{2}\times n}$ be two
sparse embedding matrices (Definition B.16) with
$t_{1}=\operatorname{poly}(r_{1}/\epsilon),t_{2}=\operatorname{poly}(r_{2}/\epsilon)$.
Then with probability at least $0.99$, $\forall X\in\mathbb{R}^{n\times s}$,

$\displaystyle(1-\epsilon)\|U\otimes V\otimes X-A\|_{F}^{2}\leq\|T_{1}U\otimes
T_{2}V\otimes X-A(T_{1},T_{2},I)\|_{F}^{2}\leq(1+\epsilon)\|U\otimes V\otimes
X-A\|_{F}^{2}.$

###### Proof.

Let $X\in\mathbb{R}^{n\times s}.$ We define $Z_{1}=(V^{\top}\odot
X^{\top})\in\mathbb{R}^{s\times n^{2}}$. We choose a sparse embedding matrix
(Definition B.16) $T_{1}\in\mathbb{R}^{t_{1}\times n}$ with
$t_{1}=\operatorname{poly}(r_{1}/\epsilon)$ rows. According to Lemma B.19 with
probability $0.999$, for all $Z\in\mathbb{R}^{s\times n^{2}}$,

$\displaystyle(1-\epsilon)\|UZ-A_{1}\|_{F}^{2}\leq\|T_{1}UZ-
T_{1}A_{1}\|_{F}^{2}\leq(1+\epsilon)\|T_{1}UZ-A_{1}\|_{F}^{2}.$

It means that

$\displaystyle(1-\epsilon)\|UZ_{1}-A_{1}\|_{F}^{2}\leq\|T_{1}UZ_{1}-T_{1}A_{1}\|_{F}^{2}\leq(1+\epsilon)\|T_{1}UZ_{1}-A_{1}\|_{F}^{2}.$

Second, we unflatten matrix $T_{1}A_{1}\in\mathbb{R}^{t_{1}\times n^{2}}$ to
obtain a tensor $A^{\prime}\in\mathbb{R}^{t_{1}\times n\times n}$. Then we
flatten $A^{\prime}$ along the second direction to obtain
$A^{\prime}_{2}\in\mathbb{R}^{n\times t_{1}n}$. We define
$Z_{2}=((T_{1}U)^{\top}\odot X^{\top})\in\mathbb{R}^{s\times t_{1}n}$. Then,
by flattening,

$\displaystyle\|V\cdot Z_{2}-A^{\prime}_{2}\|_{F}^{2}=\|T_{1}U\cdot
Z_{1}-T_{1}A_{1}\|_{F}^{2}=(1\pm\epsilon)\|U\otimes V\otimes X-A\|_{F}^{2}.$

We choose a sparse embedding matrix (Definition B.16)
$T_{2}\in\mathbb{R}^{t_{2}\times n}$ with
$t_{2}=\operatorname{poly}(r_{2}/\epsilon)$ rows. Then according to Lemma B.19
with probability $0.999$, for all $Z\in\mathbb{R}^{s\times t_{1}n}$,

$\displaystyle(1-\epsilon)\|VZ-A^{\prime}_{2}\|_{F}^{2}\leq\|T_{2}VZ-
T_{2}A^{\prime}_{2}\|_{F}^{2}\leq(1+\epsilon)\|VZ-A^{\prime}_{2}\|_{F}^{2}.$

Thus,

$\displaystyle\|T_{2}V\cdot
Z_{2}-T_{2}A^{\prime}_{2}\|_{F}^{2}=(1\pm\epsilon)^{2}\|U\otimes V\otimes
X-A\|_{F}^{2}.$

After rescaling $\epsilon$ by a constant, with probability at least $0.99$,
$\forall X\in\mathbb{R}^{n\times s}$,

$\displaystyle(1-\epsilon)\|U\otimes V\otimes X-A\|_{F}^{2}\leq\|T_{1}U\otimes
T_{2}V\otimes X-A(T_{1},T_{2},I)\|_{F}^{2}\leq(1+\epsilon)\|U\otimes V\otimes
X-A\|_{F}^{2}.$

∎

####  C.4.2 Algorithm i@

We start with a slightly unoptimized bicriteria low rank approximation
algorithm.

Algorithm 5 Frobenius Norm Bicriteria Low Rank Approximation Algorithm,
rank-$O(k^{3}/\epsilon^{3})$

1:procedure FTensorLowRankBicriteriaCubicRank($A,n,k$) $\triangleright$
Theorem C.7

2: $s_{1}\leftarrow s_{2}\leftarrow s_{3}\leftarrow O(k/\epsilon)$.

3: $t_{1}\leftarrow t_{2}\leftarrow
t_{3}\leftarrow\operatorname{poly}(k/\epsilon)$.

4: Choose $S_{i}\in\mathbb{R}^{n^{2}\times s_{i}}$ to be a Sketching matrix,
$\forall i\in[3]$. $\triangleright$ Definition B.18

5: Choose $T_{i}\in\mathbb{R}^{t_{i}\times n}$ to be a Sketching matrix,
$\forall i\in[3]$. $\triangleright$ Definition B.16

6: Compute $U\leftarrow T_{1}\cdot(A_{1}\cdot S_{1})$, $V\leftarrow
T_{2}\cdot(A_{2}\cdot S_{2})$, $W\leftarrow T_{3}\cdot(A_{3}\cdot S_{3})$.

7: Compute $C\leftarrow A(T_{1},T_{2},T_{3})$.

8:
$X\leftarrow$FTensorRegression($C,U,V,W,t_{1},s_{1},t_{2},s_{2},t_{3},s_{3}$).
$\triangleright$ Linear regression

9: return $X(A_{1}S_{1},A_{2}S_{2},A_{3}S_{3})$.

10:end procedure

######  Theorem C.7.

Given a 3rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any $k\geq
1,\epsilon\in(0,1)$, let $r=O(k^{3}/\epsilon^{3})$. There exists an algorithm
that takes $O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ time
and outputs three matrices $U\in\mathbb{R}^{n\times r},V\in\mathbb{R}^{n\times
r}$, $W\in\mathbb{R}^{n\times r}$ such that

$\displaystyle\left\|\sum_{i=1}^{r}U_{i}\otimes V_{i}\otimes
W_{i}-A\right\|_{F}^{2}\leq(1+\epsilon)\underset{\operatorname{rank}-k\leavevmode\nobreak\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2}$

holds with probability $9/10$.

###### Proof.

At the end of Theorem C.1, we need to run a polynomial system verifier. This
is why we obtain exponential in $k$ running time. Instead of running the
polynomial system verifier, we can use Lemma C.5. This reduces the running
time to be polynomial in all parameters: $n,k,1/\epsilon$. However, the output
tensor has rank $(k/\epsilon)^{3}$ (Here we mean that we do not obtain a
better decomposition than $(k/\epsilon)^{3}$ components). According to Section
B.6, for each $i$, $A_{i}S_{i}$ can be computed in
$O(\operatorname{nnz}(A))+n\operatorname{poly}(k/\epsilon)$ time. Then
$T_{i}(A_{i}S_{i})$ can be computed in $n\operatorname{poly}(k,1/\epsilon)$
time and $A(T_{1},T_{2},T_{3})$ also can be computed in
$O(\operatorname{nnz}(A))$ time. The running time for the regression is
$\operatorname{poly}(k/\epsilon)$. ∎

Now we present an optimized bicriteria algorithm.

Algorithm 6 Frobenius Norm Low Rank Approximation Algorithm,
rank-$O(k^{2}/\epsilon^{2})$

1:procedure FTensorLowRankBicriteriaQuadraticRank($A,n,k$) $\triangleright$
Theorem C.8

2: $s_{1}\leftarrow s_{2}\leftarrow O(k/\epsilon)$.

3: Choose $S_{i}\in\mathbb{R}^{n^{2}\times s_{i}}$ to be a sketching matrix,
$\forall i\in[3]$. $\triangleright$ Definition B.18

4: Compute $A_{1}\cdot S_{1}$, $A_{2}\cdot S_{2}$.

5: Form $\widehat{U}$ by using $A_{1}S_{1}$ according to Equation (9).

6: Form $\widehat{V}$ by using $A_{2}S_{2}$ according to Equation (10).

7:
$\widehat{W}\leftarrow$FTensorMultipleRegression($A,\widehat{U},\widehat{V},n,n,s_{1}s_{2}$).
$\triangleright$ Algorithm 4

8: return $\widehat{U},\widehat{V},\widehat{W}$.

9:end procedure

10:procedure FTensorLowRankBicriteriaQuadraticRank($A,n,k$) $\triangleright$
Theorem C.8

11: $s_{1}\leftarrow s_{2}\leftarrow O(k/\epsilon)$.

12: $t_{1}\leftarrow t_{2}\leftarrow\operatorname{poly}(k/\epsilon)$.

13: Choose $S_{i}\in\mathbb{R}^{n^{2}\times s_{i}}$ to be a Sketching matrix,
$\forall i\in[2]$. $\triangleright$ Definition B.18

14: Choose $T_{i}\in\mathbb{R}^{t_{i}\times n}$ to be a Sketching matrix,
$\forall i\in[2]$. $\triangleright$ Definition B.16

15: Form $\widehat{U}$ by using $A_{1}S_{1}$ according to Equation (9).

16: Form $\widehat{V}$ by using $A_{2}S_{2}$ according to Equation (10).

17: Compute $C\leftarrow A(T_{1},T_{2},I)$. $\triangleright$
$C\in\mathbb{R}^{t_{1}\times t_{2}\times n}$

18: Compute
$B\leftarrow(T_{1}\widehat{U})^{\top}\odot(T_{2}\widehat{V})^{\top}$.

19: $\widehat{W}\leftarrow\underset{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}{\arg\min}\|XB-C_{3}\|_{F}^{2}.$

20: return $\widehat{U},\widehat{V},\widehat{W}$.

21:end procedure

######  Theorem C.8.

Given a 3rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any $k\geq
1,\epsilon\in(0,1)$, let $r=O(k^{2}/\epsilon^{2})$. There exists an algorithm
that takes $O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ time
and outputs three matrices $U\in\mathbb{R}^{n\times r},V\in\mathbb{R}^{n\times
r}$, $W\in\mathbb{R}^{n\times r}$ such that

$\displaystyle\left\|\sum_{i=1}^{r}U_{i}\otimes V_{i}\otimes
W_{i}-A\right\|_{F}^{2}\leq(1+\epsilon)\underset{\operatorname{rank}-k\leavevmode\nobreak\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2}$

holds with probability $9/10$.

Note that there are two different ways to implement algorithm
FTensorLowRankBicriteriaQuadraticRank. We present the proofs for both of them
here.

Approach<EMAIL_ADDRESS>
###### Proof.

Let $\operatorname{OPT}=\underset{\operatorname{rank}-k\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2}.$ According to Theorem C.1, we know that there
exists a sketching matrix $S_{3}\in\mathbb{R}^{n^{2}\times s_{3}}$ where
$s_{3}=O(k/\epsilon)$, such that

$\displaystyle\min_{X_{1}\in\mathbb{R}^{s_{1}\times
k},X_{2}\in\mathbb{R}^{s_{2}\times k},X_{3}\in\mathbb{R}^{s_{3}\times
k}}\left\|\sum_{l=1}^{k}(A_{1}S_{1}X_{1})_{l}\otimes(A_{2}S_{2}X_{2})_{l}\otimes(A_{3}S_{3}X_{3})_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}$

Now we fix an $l$ and we have:

$\displaystyle(A_{1}S_{1}X_{1})_{l}\otimes(A_{2}S_{2}X_{2})_{l}\otimes(A_{3}S_{3}X_{3})_{l}$
$\displaystyle=$
$\displaystyle\left(\sum_{i=1}^{s_{1}}(A_{1}S_{1})_{i}(X_{1})_{i,l}\right)\otimes\left(\sum_{j=1}^{s_{2}}(A_{2}S_{2})_{j}(X_{2})_{j,l}\right)\otimes(A_{3}S_{3}X_{3})_{l}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}(A_{1}S_{1})_{i}\otimes(A_{2}S_{2})_{j}\otimes(A_{3}S_{3}X_{3})_{l}(X_{1})_{i,l}(X_{2})_{j,l}$

Thus, we have

$\displaystyle\min_{X_{1},X_{2},X_{3}}\left\|\sum_{i=1}^{s_{1}}\sum_{j=1}^{s_{2}}(A_{1}S_{1})_{i}\otimes(A_{2}S_{2})_{j}\otimes\left(\sum_{l=1}^{k}(A_{3}S_{3}X_{3})_{l}(X_{1})_{i,l}(X_{2})_{j,l}\right)-A\right\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}.$
(7)

We use matrices $A_{1}S_{1}\in\mathbb{R}^{n\times s_{1}}$ and
$A_{2}S_{2}\in\mathbb{R}^{n\times s_{2}}$ to construct a matrix
$B\in\mathbb{R}^{s_{1}s_{2}\times n^{2}}$ in the following way: each row of
$B$ is the vector corresponding to the matrix generated by the $\otimes$
product between one column vector in $A_{1}S_{1}$ and the other column vector
in $A_{2}S_{2}$, i.e.,

$\displaystyle
B^{i+(j-1)s_{1}}=\mathrm{vec}((A_{1}S_{1})_{i}\otimes(A_{2}S_{2})_{j}),\forall
i\in[s_{1}],j\in[s_{2}],$ (8)

where $(A_{1}S_{1})_{i}$ denotes the $i$-th column of $A_{1}S_{1}$ and
$(A_{2}S_{2})_{j}$ denote the $j$-th column of $A_{2}S_{2}$.

We create matrix $\widehat{U}\in\mathbb{R}^{n\times s_{1}s_{2}}$ by copying
matrix $A_{1}S_{1}$ $s_{2}$ times, i.e.,

$\displaystyle\widehat{U}=\begin{bmatrix}A_{1}S_{1}&A_{1}S_{1}&\cdots&A_{1}S_{1}\end{bmatrix}.$
(9)

We create matrix $\widehat{V}\in\mathbb{R}^{n\times s_{1}s_{2}}$ by copying
the $i$-th column of $A_{2}S_{2}$ a total of $s_{1}$ times, into columns
$(i-1)s_{1},\cdots,is_{1}$ of $\widehat{V}$, for each $i\in[s_{2}]$, i.e.,

$\displaystyle\widehat{V}=\begin{bmatrix}(A_{2}S_{2})_{1}&\cdots&(A_{2}S_{2})_{1}&(A_{2}S_{2})_{2}&\cdots&(A_{2}S_{2})_{2}&\cdots&(A_{2}S_{2})_{s_{2}}&\cdots&(A_{2}S_{2})_{s_{2}}\end{bmatrix}.$
(10)

Thus, we can use $\widehat{U}$ and $\widehat{V}$ to represent $B$,

$\displaystyle
B=(\widehat{U}^{\top}\odot\widehat{V}^{\top})\in\mathbb{R}^{s_{1}s_{2}\times
n^{2}}.$

According to Equation (7), we have:

$\displaystyle\underset{W\in\mathbb{R}^{n\times s_{1}s_{2}}}{\min}\|WB-
A_{3}\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}.$

Next, we want to find matrix $W\in\mathbb{R}^{n\times s_{1}s_{2}}$ by solving
the following optimization problem,

$\displaystyle\underset{W\in\mathbb{R}^{n\times s_{1}s_{2}}}{\min}\|WB-
A_{3}\|_{F}^{2}.$

Note that $B$ has size $s_{1}s_{2}\times n^{2}$. Naïvely writing down $B$
already requires $\Omega(n^{2})$ time. In order to achieve nearly linear time
in $n$, we cannot write down $B$. We choose
$S_{3}\in\mathbb{R}^{n_{1}n_{2}\times s_{3}}$ to be a TensorSketch (Definition
B.34). In order to solve multiple regression, we need to set
$s_{3}=O((s_{1}s_{2})^{2}+(s_{1}s_{2})/\epsilon)$. Let $\widehat{W}$ denote
the optimal solution to $\|WBS_{3}-A_{3}S_{3}\|_{F}^{2}$. Then
$\widehat{W}=(A_{3}S_{3})(BS_{3})^{\dagger}$. Since each row of $S_{3}$ has
exactly $1$ nonzero entry, $A_{3}S_{3}$ can be computed in
$O(\operatorname{nnz}(A))$ time. Since
$B=(\widehat{U}^{\top}\odot\widehat{V}^{\top})$, according to Definition B.34,
$BS_{3}$ can be computed in
$n\operatorname{poly}(s_{1}s_{2}/\epsilon)=n\operatorname{poly}(k/\epsilon)$
time. By Theorem C.4, we have

$\displaystyle\|\widehat{W}B-A_{3}\|_{F}^{2}\leq(1+\epsilon)\min_{W\in\mathbb{R}^{n\times
s_{1}s_{2}}}\|WB-A_{3}\|_{F}^{2}.$

Thus, we have

$\displaystyle\|\widehat{U}\otimes\widehat{V}\otimes\widehat{W}-A\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}.$

According to Definition B.18, $A_{1}S_{1},A_{2}S_{2}$ can be computed in
$O(\operatorname{nnz}(A)+\operatorname{poly}(k/\epsilon))$ time. Te total
running time is thus
$O(\operatorname{nnz}(A)+\operatorname{poly}(k/\epsilon))$. ∎

Approach<EMAIL_ADDRESS>
###### Proof.

Let $\operatorname{OPT}=\underset{\operatorname{rank}-k\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2}.$ Choose sketching matrices (Definition B.18)
$S_{1}\in\mathbb{R}^{n^{2}\times s_{1}}$, $S_{2}\in\mathbb{R}^{n^{2}\times
s_{2}}$, $S_{3}\in\mathbb{R}^{n^{2}\times s_{3}}$, and sketching matrices
(Definition B.16) $T_{1}\in\mathbb{R}^{t_{1}\times n}$ and
$T_{2}\in\mathbb{R}^{t_{2}\times n}$ with
$s_{1}=s_{2}=s_{3}=O(k/\epsilon),t_{1}=t_{2}=\operatorname{poly}(k/\epsilon)$.
We create matrix $\widehat{U}\in\mathbb{R}^{n\times s_{1}s_{2}}$ by copying
matrix $A_{1}S_{1}$ $s_{2}$ times, i.e.,

$\displaystyle\widehat{U}=\begin{bmatrix}A_{1}S_{1}&A_{1}S_{1}&\cdots&A_{1}S_{1}\end{bmatrix}.$

We create matrix $\widehat{V}\in\mathbb{R}^{n\times s_{1}s_{2}}$ by copying
the $i$-th column of $A_{2}S_{2}$ a total of $s_{1}$ times, into columns
$(i-1)s_{1},\cdots,is_{1}$ of $\widehat{V}$, for each $i\in[s_{2}]$, i.e.,

$\displaystyle\widehat{V}=\begin{bmatrix}(A_{2}S_{2})_{1}&\cdots&(A_{2}S_{2})_{1}&(A_{2}S_{2})_{2}&\cdots&(A_{2}S_{2})_{2}&\cdots&(A_{2}S_{2})_{s_{2}}&\cdots&(A_{2}S_{2})_{s_{2}}\end{bmatrix}.$

As we proved in Approach i@, we have

$\displaystyle\min_{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}\|\widehat{U}\otimes\widehat{V}\otimes
X-A\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}.$

Let
$B=((T_{1}\widehat{U})^{\top}\odot(T_{2}\widehat{V})^{\top})\in\mathbb{R}^{s_{1}s_{2}\times
t_{1}t_{2}},$ and flatten $A(T_{1},T_{2},I)$ along the third direction to
obtain $C_{3}\in\mathbb{R}^{n\times t_{1}t_{2}}$. Let

$\displaystyle\widehat{W}=\underset{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}{\arg\min}\|T_{1}\widehat{U}\otimes T_{2}\widehat{V}\otimes
X-A(T_{1},T_{2},I)\|_{F}^{2}=\underset{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}{\arg\min}\|XB-C_{3}\|_{F}^{2}.$

Let

$\displaystyle W^{*}=\underset{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}{\arg\min}\|\widehat{U}\otimes\widehat{V}\otimes X-A\|_{F}^{2}.$

According to Lemma C.6,

$\displaystyle\|\widehat{U}\otimes\widehat{V}\otimes\widehat{W}-A\|_{F}^{2}$
$\displaystyle\leq$
$\displaystyle\frac{1}{1-\epsilon}\|T_{1}\widehat{U}\otimes
T_{2}\widehat{V}\otimes\widehat{W}-A(T_{1},T_{2},I)\|_{F}^{2}$
$\displaystyle\leq$
$\displaystyle\frac{1}{1-\epsilon}\|T_{1}\widehat{U}\otimes
T_{2}\widehat{V}\otimes W^{*}-A(T_{1},T_{2},I)\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\frac{1+\epsilon}{1-\epsilon}\|\widehat{U}\otimes\widehat{V}\otimes
W^{*}-A\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\frac{(1+\epsilon)^{2}}{1-\epsilon}\operatorname{OPT}.$

According to Definition B.18, $A_{1}S_{1},A_{2}S_{2}$ can be computed in
$O(\operatorname{nnz}(A)+\operatorname{poly}(k/\epsilon))$ time. The total
running time is thus
$O(\operatorname{nnz}(A)+\operatorname{poly}(k/\epsilon))$. Since
$T_{1},T_{2}$ are sparse embedding matrices,
$T_{1}\widehat{U},T_{2}\widehat{V}$ can be computed in
$O(\operatorname{nnz}(A)+\operatorname{poly}(k/\epsilon))$ time. The total
running time is in $O(\operatorname{nnz}(A)+\operatorname{poly}(k/\epsilon))$.
∎

######  Theorem C.9.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any
$k\geq 1$ and any $0<\epsilon<1$, if $A_{k}$ exists then there is a randomized
algorithm running in
$\operatorname{nnz}(A)+n\cdot\operatorname{poly}(k/\epsilon)$ time which
outputs a rank-$O(k^{2}/\epsilon^{2})$ tensor $B$ for which
$\|A-B\|_{F}^{2}\leq(1+\epsilon)\|A-A_{k}\|_{F}^{2}$. If $A_{k}$ does not
exist, then the algorithm outputs a rank-$O(k^{2}/\epsilon^{2})$ tensor $B$
for which $\|A-B\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}+\gamma$, where
$\gamma$ is an arbitrarily small positive function of $n$. In both cases, the
algorithm succeeds with probability at least $9/10$.

###### Proof.

If $A_{k}$ exists, then the proof directly follows the proof of Theorem C.1
and Theorem C.8. If $A_{k}$ does not exist, then for any $\gamma>0,$ there
exist $U^{*}\in\mathbb{R}^{n\times k},V^{*}\in\mathbb{R}^{n\times
k},W^{*}\in\mathbb{R}^{n\times k}$ such that

$\displaystyle\|U^{*}\otimes V^{*}\otimes
W^{*}-A\|_{F}^{2}\leq\inf_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A-A^{\prime}\|_{F}^{2}+\frac{1}{10}\gamma.$

Then we just regard $U^{*}\otimes V^{*}\otimes W^{*}$ as the “best”
$\operatorname{rank}\leavevmode\nobreak\ k$ approximation to $A$, and follow
the same argument as in the proof of Theorem C.1 and the proof of Theorem C.8.
We can finally output a tensor $B\in\mathbb{R}^{n\times n\times n}$ with
rank-$O(k^{2}/\epsilon^{2})$ such that

$\displaystyle\|B-A\|_{F}^{2}$ $\displaystyle\leq\leavevmode\nobreak\
(1+\epsilon)\|U^{*}\otimes V^{*}\otimes W^{*}-A\|_{F}^{2}$
$\displaystyle\leq\leavevmode\nobreak\
(1+\epsilon)\left(\inf_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A-A^{\prime}\|_{F}^{2}+\frac{1}{10}\gamma\right)$
$\displaystyle\leq\leavevmode\nobreak\
(1+\epsilon)\inf_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A-A^{\prime}\|_{F}^{2}+\gamma$

where the first inequality follows by the proof of Theorem C.1 and the proof
of theorem C.8. The second inequality follows by our choice of
$U^{*},V^{*},W^{*}$. The third inequality follows since $1+\epsilon<2$ and
$\gamma>0$. ∎

####  C.4.3 $\operatorname{poly}(k)$-approximation to multiple regression

######  Lemma C.10 ((1.4) and (1.9) in [RV09]).

Let $s\geq k$. Let $U\in\mathbb{R}^{n\times k}$ denote a matrix that has
orthonormal columns, and $S\in\mathbb{R}^{s\times n}$ denote an i.i.d.
$N(0,1/s)$ Gaussian matrix. Then $SU$ is also an $s\times k$ i.i.d. Gaussian
matrix with each entry draw from $N(0,1/s)$, and furthermore, we have with
arbitrarily large constant probability,

$\displaystyle\sigma_{\max}(SU)=O(1)\text{\leavevmode\nobreak\
and\leavevmode\nobreak\ }\sigma_{\min}(SU)=\Omega(1/\sqrt{s}).$

###### Proof.

Note that
$\sqrt{s}-\sqrt{k-1}=\frac{s-k-1}{\sqrt{s}+\sqrt{k-1}}=\Omega(1/\sqrt{s})$. ∎

######  Lemma C.11.

Given matrices $A\in\mathbb{R}^{n\times k}$, $B\in\mathbb{R}^{n\times d}$, let
$S\in\mathbb{R}^{s\times n}$ denote a standard Gaussian $N(0,1)$ matrix with
$s=k$. Let $X^{*}=\underset{X\in\mathbb{R}^{k\times d}}{\min}\|AX-B\|_{F}$.
Let $X^{\prime}=\underset{X\in\mathbb{R}^{k\times d}}{\min}\|SAX-SB\|_{F}$.
Then, we have that

$\displaystyle\|AX^{\prime}-B\|_{F}\leq O(\sqrt{k})\|AX^{*}-B\|_{F},$

holds with probability at least $0.99$.

###### Proof.

Let $X^{*}\in\mathbb{R}^{k\times d}$ denote the optimal solution such that

$\displaystyle\|AX^{*}-B\|_{F}=\underset{X\in\mathbb{R}^{k\times
d}}{\min}\|AX-B\|_{F}.$

Consider a standard Gaussian matrix $S\in\mathbb{R}^{k\times n}$ scaled by
$1/\sqrt{k}$ with exactly $k$ rows. Then for any $X\in\mathbb{R}^{k\times d}$,
by the triangle inequality, we have

$\displaystyle\|SAX-SB\|_{F}\leq\|SAX-SAX^{*}\|_{F}+\|SAX^{*}-SB\|_{F},$

and

$\displaystyle\|SAX-SB\|_{F}\geq\|SAX-SAX^{*}\|_{F}-\|SAX^{*}-SB\|_{F}.$

We first show how to bound $\|SAX-SAX^{*}\|_{F}$, and then show how to bound
$\|SAX^{*}-SB\|_{F}$.

Note that Lemma C.10 implies the following result,

######  Claim C.12.

For any $X\in\mathbb{R}^{k\times d}$, with probability $0.999$, we have

$\displaystyle\frac{1}{\sqrt{k}}\|AX-AX^{*}\|_{F}\lesssim\|SAX-
SAX^{*}\|_{F}\lesssim\|AX-AX^{*}\|_{F}.$

###### Proof.

First, we can write $A=UR\in\mathbb{R}^{n\times k}$ where
$U\in\mathbb{R}^{n\times k}$ has orthonormal columns and
$R\in\mathbb{R}^{k\times k}$. It gives,

$\displaystyle\|SAX-SAX^{*}\|_{F}=\|SU(RX-RX^{*})\|_{F}.$

Second, applying Lemma C.10 to $SU\in\mathbb{R}^{s\times k}$ completes the
proof. ∎

Using Markov’s inequality, for any fixed matrix $AX^{*}-B$, choosing a
Gaussian matrix $S$, we have that

$\displaystyle\|SAX^{*}-SB\|_{F}^{2}=O(\|AX^{*}-B\|_{F}^{2})$

holds with probability at least $0.999$. This is equivalent to

$\displaystyle\|SAX^{*}-SB\|_{F}=O(\|AX^{*}-B\|_{F}),$ (11)

holding with probability at least $0.999$.

Let $X^{\prime}=\underset{X\in\mathbb{R}^{k\times d}}{\arg\min}\|SAX-
SB\|_{F}$. Putting it all together, we have

$\displaystyle\leavevmode\nobreak\ \|AX^{\prime}-B\|_{F}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\|AX^{\prime}-AX^{*}\|_{F}+\|AX^{*}-B\|_{F}$ by triangle inequality
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
O(\sqrt{k})\|SAX^{\prime}-SAX^{*}\|_{F}+\|AX^{*}-B\|_{F}$ by Claim C.12
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
O(\sqrt{k})\|SAX^{\prime}-SB\|_{F}+O(\sqrt{k})\|SAX^{*}-SB\|_{F}+\|AX^{*}-B\|_{F}$
by triangle inequality $\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
O(\sqrt{k})\|SAX^{*}-SB\|_{F}+O(\sqrt{k})\|SAX^{*}-SB\|_{F}+\|AX^{*}-B\|_{F}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
definition\leavevmode\nobreak\ of\leavevmode\nobreak\ }X^{\prime}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
O(\sqrt{k})\|AX^{*}-B\|_{F}.$ by Equation (11)

∎

####  C.4.4 Algorithm ii@

######  Theorem C.13.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any
$k\geq 1$, let $r=k^{2}$. There exists an algorithm which takes
$O(\operatorname{nnz}(A)k)+n\operatorname{poly}(k)$ time and outputs three
matrices $U,V,W\in\mathbb{R}^{n\times r}$ such that,

$\displaystyle\left\|\sum_{i=1}^{r}U_{i}\otimes V_{i}\otimes
W_{i}-A\right\|_{F}\leq\operatorname{poly}(k)\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}$

holds with probability $9/10$.

###### Proof.

Let $\operatorname{OPT}=\underset{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}{\min}\|A^{\prime}-A\|_{F},$ we fix $V^{*}\in\mathbb{R}^{n\times
k},W^{*}\in\mathbb{R}^{n\times k}$ to be the optimal solution of the original
problem. We use $Z_{1}=(V^{*\top}\odot W^{*\top})\in\mathbb{R}^{k\times
n^{2}}$ to denote the matrix where the $i$-th row is the vectorization of
$V_{i}^{*}\otimes W_{i}^{*}$. Let $A_{1}\in\mathbb{R}^{n\times n^{2}}$ denote
the matrix obtained by flattening tensor $A\in\mathbb{R}^{n\times n\times n}$
along the first direction. Then, we have

$\displaystyle\min_{U}\|UZ_{1}-A_{1}\|_{F}\leq\operatorname{OPT}.$

Choosing an $N(0,1/k)$ Gaussian sketching matrix
$S_{1}\in\mathbb{R}^{n^{2}\times s_{1}}$ with $s_{1}=k$, we can obtain the
smaller problem,

$\displaystyle\min_{U\in\mathbb{R}^{n\times
k}}\|UZ_{1}S_{1}-A_{1}S_{1}\|_{F}.$

Define $\widehat{U}=A_{1}S_{1}(Z_{1}S_{1})^{\dagger}$. Define
$\alpha=O(\sqrt{k})$. By Lemma C.11, we have

$\displaystyle\|\widehat{U}Z_{1}-A_{1}\|_{F}\leq\alpha\operatorname{OPT}.$

Second, we fix $\widehat{U}$ and $W^{*}$. Define $Z_{2},A_{2}$ similarly as
above. Choosing an $N(0,1/k)$ Gaussian sketching matrix
$S_{2}\in\mathbb{R}^{n^{2}\times s_{2}}$ with $s_{2}=k$, we can obtain another
smaller problem,

$\displaystyle\min_{V\in\mathbb{R}^{n\times
k}}\|VZ_{2}S_{2}-A_{2}S_{2}\|_{F}.$

Define $\widehat{V}=A_{2}S_{2}(Z_{2}S_{2})^{\dagger}$. By Lemma C.11 again, we
have

$\displaystyle\|\widehat{V}Z_{2}-A_{2}\|_{F}\leq\alpha^{2}\operatorname{OPT}.$

Thus, we now have

$\displaystyle\min_{X_{1},X_{2},W}\|A_{1}S_{1}X_{1}\otimes
A_{2}S_{2}X_{2}\otimes W-A\|_{F}\leq\alpha^{2}\operatorname{OPT}$

We use a similar idea as in the proof of Theorem C.8. We create matrix
$\widetilde{U}\in\mathbb{R}^{n\times s_{1}s_{2}}$ by copying matrix
$A_{1}S_{1}$ $s_{2}$ times, i.e.,

$\displaystyle\widetilde{U}=\begin{bmatrix}A_{1}S_{1}&A_{1}S_{1}&\cdots&A_{1}S_{1}\end{bmatrix}.$

We create matrix $\widetilde{V}\in\mathbb{R}^{n\times s_{1}s_{2}}$ by copying
the $i$-th column of $A_{2}S_{2}$ a total of $s_{1}$ times, into columns
$(i-1)s_{1},\cdots,is_{1}$ of $\widetilde{V}$, for each $i\in[s_{2}]$, i.e.,

$\displaystyle\widetilde{V}=\begin{bmatrix}(A_{2}S_{2})_{1}&\cdots&(A_{2}S_{2})_{1}&(A_{2}S_{2})_{2}&\cdots&(A_{2}S_{2})_{2}&\cdots&(A_{2}S_{2})_{s_{2}}&\cdots&(A_{2}S_{2})_{s_{2}}\end{bmatrix}.$

We have

$\displaystyle\min_{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}\|\widetilde{U}\otimes\widetilde{V}\otimes
X-A\|_{F}\leq\alpha^{2}\operatorname{OPT}.$

Choose $T_{i}\in\mathbb{R}^{t_{i}\times n}$ to be a sparse embedding matrix
(Definition B.16) with $t_{i}=\operatorname{poly}(k/\epsilon)$, for each
$i\in[2]$. By applying Lemma C.6, we have, if $W^{\prime}$ satisfies,

$\displaystyle\|T_{1}\widetilde{U}\otimes T_{2}\widetilde{V}\otimes
W^{\prime}-A(T_{1},T_{2},I)\|_{F}=\min_{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}\|T_{1}\widetilde{U}\otimes T_{2}\widetilde{V}\otimes
X-A(T_{1},T_{2},I)\|_{F}$

then,

$\displaystyle\|\widetilde{U}\otimes\widetilde{V}\otimes
W^{\prime}-A\|_{F}\leq(1+\epsilon)\min_{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}\|\widetilde{U}\otimes\widetilde{V}\otimes
X-A\|_{F}\leq(1+\epsilon)\alpha^{2}\operatorname{OPT}.$

Thus, we only need to solve

$\displaystyle\min_{X\in\mathbb{R}^{n\times
s_{1}s_{2}}}\|T_{1}\widetilde{U}\otimes T_{2}\widetilde{V}\otimes
X-A(T_{1},T_{2},I)\|_{F}.$

which is similar to the proof of Theorem C.8. Therefore, we complete the proof
of correctness. For the running time, $A_{1}S_{1},A_{2}S_{2}$ can be computed
in $O(\operatorname{nnz}(A)k)$ time, $T_{1}\widetilde{U},T_{2}\widetilde{V}$
can be computed in $n\operatorname{poly}(k)$ time. The final regression
problem can be computed in $n\operatorname{poly}(k)$ running time. ∎

###  C.5 Generalized matrix row subset selection

Note that in this section, the notation $\Pi_{C,k}^{\xi}$ is given in
Definition B.5.

Algorithm 7 Generalized Matrix Row Subset Selection: Constructing $R$ with
$r=O(k+k/\epsilon)$ Rows and a rank-$k$ $U\in\mathbb{R}^{k\times r}$

1:procedure GeneralizedMatrixRowSubsetSelection($A,C,n,m,k,\epsilon$)
$\triangleright$ Theorem C.14

2: $Y,\Phi,\Delta\leftarrow$ ApproxSubspaceSVD($A,C,k$). $\triangleright$
Claim C.16 and Lemma 3.12 in [BW14]

3: $B\leftarrow Y\Delta$.

4: $Z_{2},D\leftarrow\textsc{QR}(B)$. $\triangleright$
$Z_{2}\in\mathbb{R}^{m\times k}$, $Z_{2}^{\top}Z_{2}=I_{k}$,
$D\in\mathbb{R}^{k\times k}$

5: $h_{2}\leftarrow 8k\ln(20k)$.

6: $\Omega_{2},D_{2}\leftarrow\textsc{RandSampling}(Z_{2},h_{2},1)$
$\triangleright$ Definition 3.6 in [BW14]

7: $M_{2}\leftarrow Z_{2}^{\top}\Omega_{2}D_{2}\in\mathbb{R}^{k\times h_{2}}$.

8: $U_{M_{2}},\Sigma_{M_{2}},V_{M_{2}}^{\top}\leftarrow\textsc{SVD}(M_{2})$.
$\triangleright$ $\operatorname{rank}(M_{2})=k$ and
$V_{M_{2}}\in\mathbb{R}^{h_{2}\times k}$

9: $r_{1}\leftarrow 4k$.

10: $S_{2}\leftarrow$
BSSSamplingSparse($V_{M_{2}},((A^{\top}-A^{\top}Z_{2}Z_{2}^{\top})\Omega_{2}D_{2})^{\top},r_{1},0.5$)
$\triangleright$ Lemma 4.3 in [BW14]

11:
$R_{1}\leftarrow(A^{\top}\Omega_{2}D_{2}S_{2})^{\top}\in\mathbb{R}^{r_{1}\times
n}$ containing rescaled rows from $A$.

12: $r_{2}\leftarrow 4820k/\epsilon$.

13: $R_{2}\leftarrow$ AdaptiveRowsSparse($A,Z_{2},R_{1},r_{2}$)
$\triangleright$ Lemma 4.5 in [BW14]

14: $R\leftarrow[R_{1}^{\top},R_{2}^{\top}]^{\top}$. $\triangleright$
$R\in\mathbb{R}^{(r_{1}+r_{2})\times n}$ containing $r=4k+4820k/\epsilon$
rescaled rows of $A$.

15: Choose $W\in\mathbb{R}^{\xi\times m}$ to be a randomly chosen sparse
subspace embedding with $\xi=\Omega(k^{2}\epsilon^{-2})$.

16: $U\leftarrow\Phi^{-1}\Delta D^{-1}(WC\Phi^{-1}\Delta
D^{-1})^{\dagger}WAR^{\dagger}=\Phi^{-1}\Delta\Delta^{\top}(WC)^{\dagger}WAR^{\dagger}$.

17: return $R$, $U$.

18:end procedure

######  Theorem C.14.

Given matrices $A\in\mathbb{R}^{n\times m}$ and $C\in\mathbb{R}^{n\times k}$,
there exists an algorithm which takes $O(\operatorname{nnz}(A)\log
n)+(m+n)\operatorname{poly}(k,1/\epsilon)$ time and outputs a diagonal matrix
$D\in\mathbb{R}^{n\times n}$ with $d=O(k/\epsilon)$ nonzeros (or equivalently
a matrix $R$ that contains $d=O(k/\epsilon)$ rescaled rows of $A$) and a
matrix $U\in\mathbb{R}^{k\times d}$ such that

$\displaystyle\|CUDA-A\|_{F}^{2}\leq(1+\epsilon)\min_{X\in\mathbb{R}^{k\times
m}}\|CX-A\|_{F}^{2}$

holds with probability $.99$.

###### Proof.

This follows by combining Lemma C.17 and C.18. Let $U,R$ denote the output of
procedure GeneralizedMatrixRowSubsetSelection,

$\displaystyle\|A-CUR\|_{F}^{2}\leq$ $\displaystyle\leavevmode\nobreak\
(1+\epsilon)\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
(1+\epsilon)(1+60\epsilon)\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
(1+130\epsilon)\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}.$

Because $R$ is a subset of rows of $A$ and $R$ has size $O(k/\epsilon)\times
m$, there must exist a diagonal matrix $D\in\mathbb{R}^{n\times n}$ with
$O(k/\epsilon)$ nonzeros such that $R=DA$. This completes the proof. ∎

######  Corollary C.15 (A slightly different version of Theorem C.14, faster
running time, and small input matrix).

Given matrices $A\in\mathbb{R}^{n\times m}$ and $C\in\mathbb{R}^{n\times k}$,
if $\min(m,n)=\operatorname{poly}(k,1/\epsilon)$, then there exists an
algorithm which takes
$O(\operatorname{nnz}(A))+(m+n)\operatorname{poly}(k,1/\epsilon)$ time and
outputs a diagonal matrix $D\in\mathbb{R}^{n\times n}$ with $d=O(k/\epsilon)$
nonzeros (or equivalently a matrix $R$ that contains $d=O(k/\epsilon)$
rescaled rows of $A$) and a matrix $U\in\mathbb{R}^{k\times d}$ such that

$\displaystyle\|CUDA-A\|_{F}^{2}\leq(1+\epsilon)\min_{X\in\mathbb{R}^{k\times
m}}\|CX-A\|_{F}^{2}$

holds with probability $.99$.

###### Proof.

The $\log n$ factor comes from the adaptive sampling where we need to choose a
Gaussian matrix with $O(\log n)$ rows and compute $SA$. If $A$ has
$\operatorname{poly}(k,1/\epsilon)$ columns, it is sufficient to choose $S$ to
be a CountSketch matrix with $\operatorname{poly}(k,1/\epsilon)$ rows. Then,
we do not need a $\log n$ factor in the running time. If $S$ has
$\operatorname{poly}(k,1/\epsilon)$ rows, then we no longer need the matrix
$S$. ∎

######  Claim C.16.

Given matrices $A\in\mathbb{R}^{m\times n}$ and $C\in\mathbb{R}^{m\times c}$,
let $Y\in\mathbb{R}^{m\times c},\Phi\in\mathbb{R}^{c\times c}$ and
$\Delta\in\mathbb{R}^{c\times k}$ denote the output of procedure
$\textsc{ApproxSubspaceSVD}(A,C,k,\epsilon)$. Then with probability $.99$, we
have,

$\displaystyle\|A-Y\Delta\Delta^{\top}Y^{\top}A\|_{F}^{2}\leq(1+30\epsilon)\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}.$

###### Proof.

This follows by Lemma 3.12 in [BW14]. ∎

######  Lemma C.17.

The matrices $R$ and $Z_{2}$ in procedure GeneralizedMatrixRowSubsetSelection
(Algorithm 7) satisfy with probability at least $0.17-2/n$,

$\displaystyle\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}\leq\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}+60\epsilon\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}.$

###### Proof.

We can show,

$\displaystyle\leavevmode\nobreak\
\|A-Z_{2}Z_{2}^{\top}A\|_{F}^{2}+\frac{30\epsilon}{4820}\|A-AR_{1}^{\dagger}R_{1}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|A-BB^{\dagger}A\|_{F}^{2}+\frac{30\epsilon}{4820}\|A-AR_{1}^{\dagger}R_{1}\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\|A-BB^{\dagger}A\|_{F}^{2}+30\epsilon\|A-A_{k}\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\|A-Y\Delta\Delta^{\top}YA\|_{F}^{2}+30\epsilon\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
(1+30\epsilon)\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2}+30\epsilon\|A-\Pi_{C,k}^{F}(A)\|_{F}^{2},$

where the first step follows by the fact that
$Z_{2}Z_{2}^{\top}=Z_{2}DD^{-1}Z_{2}^{\top}=(Z_{2}D)(Z_{2}D)^{\dagger}=BB^{\dagger}$,
the second step follows by $\|A-AR_{1}^{\dagger}R_{1}\|_{F}^{2}\leq
4820\|A-A_{k}\|_{F}^{2}$, the third step follows by $B=Y\Delta$ and
$B^{\dagger}=(Y\Delta)^{\dagger}=\Delta^{\dagger}Y^{\dagger}=\Delta^{\top}Y^{\top}$,
and the last step follows by Claim C.16. ∎

######  Lemma C.18.

The matrices $C,U$ and $R$ in procedure GeneralizedMatrixRowSubsetSelection
(Algorithm 7) satisfy that

$\displaystyle\|A-CUR\|_{F}^{2}\leq(1+\epsilon)\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}$

with probability at least $.99$.

###### Proof.

Let $U_{R},\Sigma_{R},V_{R}$ denote the SVD of $R$. Then
$V_{R}V_{R}^{\top}=R^{\dagger}R$.

We define $Y^{*}$ to be the optimal solution of

$\displaystyle\min_{X\in\mathbb{R}^{k\times
r}}\|WAV_{R}V_{R}^{\top}-WC\Phi^{-1}\Delta D^{-1}YR\|_{F}^{2}.$

We define $\widehat{X}^{*}$ to be $Y^{*}R\in\mathbb{R}^{k\times n}$, which is
also equivalent to defining $\widehat{X}^{*}$ to be the optimal solution of

$\displaystyle\min_{X\in\mathbb{R}^{k\times
n}}\|WAV_{R}V_{R}^{\top}-WC\Phi^{-1}\Delta D^{-1}X\|_{F}^{2}.$

Furthermore, it implies $\widehat{X}^{*}=(WC\Phi^{-1}\Delta
D^{-1})^{\dagger}WAV_{R}V_{R}^{\dagger}$.

We also define $X^{*}$ to be the optimal solution of

$\displaystyle\min_{X\in\mathbb{R}^{k\times
n}}\|AV_{R}V_{R}^{\dagger}-C\Phi^{-1}\Delta D^{-1}X\|_{F}^{2},$

which implies that,

$\displaystyle X^{*}=(C\Phi^{-1}\Delta
D^{-1})^{\dagger}AV_{R}V_{R}^{\top}=Z_{2}^{\top}AV_{R}V_{R}^{\top}.$

Now, we start to prove an upper bound on $\|A-CUR\|_{F}^{2}$,

$\displaystyle\|A-CUR\|_{F}^{2}=$ $\displaystyle\leavevmode\nobreak\
\|A-C\Phi^{-1}\Delta D^{-1}Y^{*}R\|_{F}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
definition\leavevmode\nobreak\ of\leavevmode\nobreak\ }U$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|A-C\Phi^{-1}\Delta
D^{-1}\widehat{X}^{*}\|_{F}^{2}$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }\widehat{X}^{*}=Y^{*}R$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}\widehat{X}^{*}+A-AV_{R}V_{R}^{\top}\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\underbrace{\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}\widehat{X}^{*}\|_{F}^{2}}_{\alpha}+\underbrace{\|A-AV_{R}V_{R}^{\top}\|_{F}^{2}}_{\beta},$
(12)

where the last step follows by $\widehat{X}^{*}=MV_{R}^{\top}$,
$A-AV_{R}V_{R}^{\top}=A(I-V_{R}V_{R}^{\top})$ and the Pythagorean theorem. We
show how to upper bound the term $\alpha$,

$\displaystyle\alpha\leq$ $\displaystyle\leavevmode\nobreak\
(1+\epsilon)\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta D^{-1}X^{*}\|_{F}^{2}$ by
Lemma C.19 $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\epsilon\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}X^{*}\|_{F}^{2}+\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}X^{*}\|_{F}^{2}$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\epsilon\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}X^{*}\|_{F}^{2}+\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}(Z_{2}^{\top}AR^{\dagger}R)\|_{F}^{2}.$ (13)

By the Pythagorean theorem and the definition of $Z_{2}$ (which means
$Z_{2}=C\Phi^{-1}\Delta D^{-1}$), we have,

$\displaystyle\leavevmode\nobreak\ \|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}(Z_{2}^{\top}AR^{\dagger}R)\|_{F}^{2}+\beta$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}(Z_{2}^{\top}AR^{\dagger}R)\|_{F}^{2}+\|A-AV_{R}V_{R}^{\top}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\ \|A-C\Phi^{-1}\Delta
D^{-1}(Z_{2}^{\top}AR^{\dagger}R)\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}.$ (14)

Combining Equations (C.5), (C.5) and (C.5) together, we obtain,

$\displaystyle\|A-CUR\|_{F}^{2}\leq\epsilon\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}X^{*}\|_{F}^{2}+\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}.$

We want to show $\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}X^{*}\|_{F}^{2}\leq\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}$,

$\displaystyle\leavevmode\nobreak\ \|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}X^{*}\|_{F}^{2}$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|AV_{R}V_{R}^{\top}-C\Phi^{-1}\Delta
D^{-1}Z_{2}^{\top}AV_{R}V_{R}^{\top}\|_{F}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}X^{*}=Z_{2}^{\top}AV_{R}V_{R}^{\top}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \|A-C\Phi^{-1}\Delta
D^{-1}Z_{2}^{\top}A\|_{F}^{2}$ by properties of projections
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\ \|A-C\Phi^{-1}\Delta
D^{-1}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}$ by properties of projections
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|A-Z_{2}Z_{2}^{\top}AR^{\dagger}R\|_{F}^{2}.$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}Z_{2}=C\Phi^{-1}\Delta D^{-1}$

This completes the proof. ∎

######  Lemma C.19 ([CW13]).

Let $A\in\mathbb{R}^{n\times d}$ have rank $\rho$ and $B\in\mathbb{R}^{n\times
r}$. Let $W\in\mathbb{R}^{r\times n}$ be a randomly chosen sparse subspace
embedding with $r=\Omega(\rho^{2}\epsilon^{-2})$. Let
$\widehat{X}^{*}=\underset{X\in\mathbb{R}^{d\times r}}{\arg\min}\|WAX-
WB\|_{F}^{2}$ and let $X^{*}=\underset{X\in\mathbb{R}^{d\times
r}}{\arg\min}\|AX-B\|_{F}^{2}$. Then with probability at least $.99$,

$\displaystyle\|A\widetilde{X}^{*}-B\|_{F}^{2}\leq(1+\epsilon)\|AX^{*}-B\|_{F}^{2}.$

###  C.6 Column, row, and tube subset selection, $(1+\epsilon)$-approximation

Algorithm 8 Frobenius Norm Tensor Column, Row and Tube Subset Selection,
Polynomial Time

1:procedure FCRTSelection($A,n,k,\epsilon$) $\triangleright$ Theorem C.20

2: $s_{1}\leftarrow s_{2}\leftarrow O(k/\epsilon)$.

3: Choose a Gaussian matrix $S_{1}$ with $s_{1}$ columns. $\triangleright$
Definition B.18

4: Choose a Gaussian matrix $S_{2}$ with $s_{2}$ columns. $\triangleright$
Definition B.18

5: Form matrix $Z_{3}^{\prime}$ by setting the $(i,j)$-th row to be the
vectorization of $(A_{1}S_{1})_{i}\otimes(A_{2}S_{2})_{j}$.

6: $D_{3}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{3}^{\top}$,
$(Z_{3}^{\prime})^{\top}$,$n^{2}$,$n$,$s_{1}s_{2}$,$\epsilon$).
$\triangleright$ Algorithm 7

7: Let $d_{3}$ denote the number of nonzero entries in $D_{3}$.
$\triangleright$ $d_{3}=O(s_{1}s_{2}/\epsilon)$

8: Form matrix $Z_{2}^{\prime}$ by setting the $(i,j)$-th row to be the
vectorization of $(A_{1}S_{1})_{i}\otimes(A_{3}S_{3}^{\prime})_{j}$.

9: $D_{2}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{2}^{\top}$,
$(Z_{2}^{\prime})^{\top}$,$n^{2}$,$n$,$s_{1}d_{3}$,$\epsilon$).

10: Let $d_{2}$ denote the number of nonzero entries in $D_{2}$.
$\triangleright$ $d_{2}=O(s_{1}d_{3}/\epsilon)$

11: Form matrix $Z_{1}^{\prime}$ by setting the $(i,j)$-th row to be the
vectorization of $(A_{2}D_{2})_{i}\otimes(A_{3}D_{3})_{j}$.

12: $D_{1}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{1}^{\top}$,
$(Z_{1}^{\prime})^{\top}$,$n^{2}$,$n$,$d_{2}d_{3}$,$\epsilon$).

13: Let $d_{1}$ denote the number of nonzero entries in $D_{1}$.
$\triangleright$ $d_{1}=O(d_{2}d_{3}/\epsilon)$

14: $C\leftarrow A_{1}D_{1}$, $R\leftarrow A_{2}D_{2}$ and $T\leftarrow
A_{3}D_{3}$.

15: return $C$, $R$ and $T$.

16:end procedure

We provide two bicriteria CURT results in this Section. We first present a
warm-up result. That result (Theorem C.20) does not output tensor $U$ and only
guarantees that there is a
$\operatorname{rank}$-$\operatorname{poly}(k/\epsilon)$ tensor $U$. Then we
show the second result (Theorem C.21), our second result is able to output
tensor $U$. The $U$ has rank $\operatorname{poly}(k/\epsilon)$, but not $k$.

######  Theorem C.20.

Given a 3rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any $k\geq
1$, there exists an algorithm which takes
$O(\operatorname{nnz}(A))+n\operatorname{poly}(k,1/\epsilon)$ time and outputs
three matrices: $C\in\mathbb{R}^{n\times c}$, a subset of columns of $A$,
$R\in\mathbb{R}^{n\times r}$ a subset of rows of $A$, and
$T\in\mathbb{R}^{n\times t}$, a subset of tubes of $A$ where
$c=r=t=\operatorname{poly}(k,1/\epsilon)$, and there exists a tensor
$U\in\mathbb{R}^{c\times r\times t}$ such that

$\displaystyle\|(((U\cdot T^{\top})^{\top}\cdot R^{\top})^{\top}\cdot
C^{\top})^{\top}-A\|_{F}^{2}\leq(1+\epsilon)\underset{\operatorname{rank}-k\leavevmode\nobreak\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2},$

or equivalently,

$\displaystyle\left\|\sum_{i=1}^{c}\sum_{j=1}^{r}\sum_{l=1}^{t}U_{i,j,l}\cdot
C_{i}\otimes R_{j}\otimes
T_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\underset{\operatorname{rank}-k\leavevmode\nobreak\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2}$

holds with probability $9/10$.

###### Proof.

We mainly analyze Algorithm 8 and it is easy to extend to Algorithm 9.

We fix $V^{*}\in\mathbb{R}^{n\times k}$ and $W^{*}\in\mathbb{R}^{n\times k}$.
We define $Z_{1}\in\mathbb{R}^{k\times n^{2}}$ where the $i$-th row of $Z_{1}$
is the vector $V_{i}\otimes W_{i}$. Choose sketching (Gaussian) matrix
$S_{1}\in\mathbb{R}^{n^{2}\times s_{1}}$ (Definition B.18), and let
$\widehat{U}=A_{1}S_{1}(Z_{1}S_{1})^{\dagger}\in\mathbb{R}^{n\times k}$.
Following a similar argument as in the previous theorem, we have

$\displaystyle\|\widehat{U}Z_{1}-A_{1}\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}.$

We fix $\widehat{U}$ and $W^{*}$. We define $Z_{2}\in\mathbb{R}^{k\times
n^{2}}$ where the $i$-th row of $Z_{2}$ is the vector $\widehat{U}_{i}\otimes
W^{*}_{i}$. Choose sketching (Gaussian) matrix
$S_{2}\in\mathbb{R}^{n^{2}\times s_{2}}$ (Definition B.18), and let
$\widehat{V}=A_{2}S_{2}(Z_{2}S_{2})^{\dagger}\in\mathbb{R}^{n\times k}$.
Following a similar argument as in the previous theorem, we have

$\displaystyle\|\widehat{V}Z_{2}-A_{2}\|_{F}^{2}\leq(1+\epsilon)^{2}\operatorname{OPT}.$

We fix $\widehat{U}$ and $\widehat{V}$. Note that
$\widehat{U}=A_{1}S_{1}(Z_{1}S_{1})^{\dagger}$ and
$\widehat{V}=A_{2}S_{2}(Z_{2}S_{2})^{\dagger}$. We define
$Z_{3}\in\mathbb{R}^{k\times n^{2}}$ such that the $i$-th row of $Z_{3}$ is
the vector $\widehat{U}_{i}\otimes\widehat{V}_{i}$. Let $z_{3}=s_{1}\cdot
s_{2}$. We define $Z^{\prime}_{3}\in\mathbb{R}^{z_{3}\times n^{2}}$ such that,
$\forall i\in[s_{1}],\forall j\in[s_{2}]$, the $i+(j-1)s_{1}$-th row of
$Z^{\prime}_{3}$ is the vector $(A_{1}S_{1})_{i}\otimes(A_{2}S_{2})_{j}$. We
consider the following objective function,

$\displaystyle\min_{W\in\mathbb{R}^{n\times k},X\in\mathbb{R}^{k\times
z_{3}}}\|WXZ_{3}^{\prime}-A_{3}\|_{F}^{2}\leq\min_{W\in\mathbb{R}^{n\times
k}}\|WZ_{3}-A_{3}\|_{F}^{2}\leq(1+\epsilon)^{2}\operatorname{OPT}.$

Using Theorem C.14, we can find a diagonal matrix
$D_{3}\in\mathbb{R}^{n^{2}\times n^{2}}$ with
$d_{3}=O(z_{3}/\epsilon)=O(k^{2}/\epsilon^{3})$ nonzero entries such that

$\displaystyle\min_{X\in\mathbb{R}^{d_{3}\times
z_{3}}}\|A_{3}D_{3}XZ_{3}^{\prime}-A_{3}\|_{F}^{2}\leq(1+\epsilon)^{3}\operatorname{OPT}.$

In the following, we abuse notation and let $A_{3}D_{3}\in\mathbb{R}^{n\times
d_{3}}$ by deleting zero columns. Let $W^{\prime}$ denote
$A_{3}D_{3}\in\mathbb{R}^{n\times d_{3}}$. Then,

$\displaystyle\min_{X\in\mathbb{R}^{d_{3}\times
z_{3}}}\|W^{\prime}XZ_{3}^{\prime}-A_{3}\|_{F}^{2}\leq(1+\epsilon)^{3}\operatorname{OPT}.$

We fix $\widehat{U}$ and $W^{\prime}$. Let $z_{2}=s_{1}\cdot d_{3}$. We define
$Z^{\prime}_{2}\in\mathbb{R}^{z_{2}\times n^{2}}$ such that, $\forall
i\in[s_{1}],\forall j\in[d_{3}]$, the $i+(j-1)s_{1}$-th row of
$Z^{\prime}_{2}$ is the vector $(A_{1}S_{1})_{i}\otimes(A_{3}D_{3})_{j}$.

Using Theorem C.14, we can find a diagonal matrix
$D_{2}\in\mathbb{R}^{n^{2}\times n^{2}}$ with
$d_{2}=O(z_{2}/\epsilon)=O(s_{1}d_{3}/\epsilon)=O(k^{3}/\epsilon^{5})$ nonzero
entries such that

$\displaystyle\min_{X\in\mathbb{R}^{d_{2}\times
z_{2}}}\|A_{2}D_{2}XZ_{2}^{\prime}-A_{2}\|_{F}^{2}\leq(1+\epsilon)^{4}\operatorname{OPT}.$

Let $V^{\prime}$ denote $A_{2}D_{2}$. Then,

$\displaystyle\min_{X\in\mathbb{R}^{d_{2}\times
z_{2}}}\|V^{\prime}XZ^{\prime}_{2}-A_{2}\|_{F}^{2}\leq(1+\epsilon)^{4}\operatorname{OPT}.$

We fix $V^{\prime}$ and $W^{\prime}$. Let $z_{1}=d_{2}\cdot d_{3}$. We define
$Z^{\prime}_{1}\in\mathbb{R}^{z_{1}\times n^{2}}$ such that, $\forall
i\in[d_{2}],\forall j\in[d_{3}]$, the $i+(j-1)s_{1}$-th row of
$Z^{\prime}_{1}$ is the vector $(A_{2}D_{2})_{i}\otimes(A_{3}D_{3})_{j}$.

Using Theorem C.14, we can find a diagonal matrix
$D_{1}\in\mathbb{R}^{n^{2}\times n^{2}}$ with
$d_{1}=O(z_{1}/\epsilon)=O(d_{2}d_{3}/\epsilon)=O(k^{5}/\epsilon^{9})$ nonzero
entries such that

$\displaystyle\min_{X\in\mathbb{R}^{d_{1}\times
z_{1}}}\|A_{1}D_{1}XZ_{1}^{\prime}-A_{1}\|_{F}^{2}\leq(1+\epsilon)^{5}\operatorname{OPT}.$

Let $U^{\prime}$ denote $A_{1}D_{1}$. Then,

$\displaystyle\min_{X\in\mathbb{R}^{d_{1}\times
z_{1}}}\|U^{\prime}XZ^{\prime}_{1}-A_{1}\|_{F}^{2}\leq(1+\epsilon)^{5}\operatorname{OPT}.$

Putting $U^{\prime},V^{\prime},W^{\prime}$ all together, we complete the
proof.

All the above analysis gives the running time $O(\operatorname{nnz}(A))\log
n+n^{2}\operatorname{poly}(\log n,k,1/\epsilon)$. To improve the running time,
we need to use Algorithm 9, the similar analysis will go through, the running
time will be improved to
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$, but the sample
complexity of $c,r,k$ will be slightly worse ($\operatorname{poly}\log$
factors). ∎

Algorithm 9 Frobenius Norm Tensor Column, Row and Tube Subset Selection, Input
Sparsity Time

1:procedure FCRTSelection($A,n,k,\epsilon$) $\triangleright$ Theorem C.20

2: $s_{1}\leftarrow s_{2}\leftarrow O(k/\epsilon)$.

3: $\epsilon_{0}\leftarrow 0.001$.

4: Choose a Gaussian matrix $S_{1}$ with $s_{1}$ columns. $\triangleright$
Definition B.18

5: Choose a Gaussian matrix $S_{2}$ with $s_{2}$ columns. $\triangleright$
Definition B.18

6: Form matrix $B_{1}$ by setting $(i,j)$-th column to be $(A_{1}S_{1})_{i}$.

7: Form matrix $B_{2}$ by setting $(i,j)$-th column to be $(A_{2}S_{2})_{j}$.
$\triangleright$ $Z_{3}^{\prime}=B_{1}^{\top}\odot B_{2}^{\top}$

8: $d_{3}\leftarrow O(s_{1}s_{2}\log(s_{1}s_{2})+(s_{1}s_{2}/\epsilon))$.

9:
$D_{3}\leftarrow$FastTensorLeverageScoreGeneralOrder($B_{1}^{\top},B_{2}^{\top},n,n,s_{1}s_{2},\epsilon_{0},d_{1}$).
$\triangleright$ Algorithm 15

10: Form matrix $B_{1}$ by setting $(i,j)$-th column to be $(A_{1}S_{1})_{i}$.

11: Form matrix $B_{3}$ by setting $(i,j)$-th column to be $(A_{3}D_{3})_{j}$.
$\triangleright$ $Z_{2}^{\prime}=B_{1}^{\top}\odot B_{3}^{\top}$

12: $d_{2}\leftarrow O(s_{1}d_{3}\log(s_{1}d_{3})+(s_{1}d_{3}/\epsilon))$.

13:
$D_{2}\leftarrow$FastTensorLeverageScoreGeneralOrder($B_{1}^{\top},B_{3}^{\top},n,n,s_{1}d_{3},\epsilon_{0},d_{2}$).

14: Form matrix $B_{2}$ by setting $(i,j)$-th column to be $(A_{2}D_{2})_{i}$.

15: Form matrix $B_{3}$ by setting $(i,j)$-th column to be $(A_{3}D_{3})_{j}$.
$\triangleright$ $Z_{1}^{\prime}=B_{2}^{\top}\odot B_{3}^{\top}$

16: $d_{1}\leftarrow O(d_{2}d_{3}\log(d_{2}d_{3})+(d_{2}d_{3}/\epsilon))$.

17:
$D_{1}\leftarrow$FastTensorLeverageScoreGeneralOrder($B_{2}^{\top},B_{3}^{\top},n,n,d_{2}d_{3},\epsilon_{0},d_{1}$).

18: $C\leftarrow A_{1}D_{1}$, $R\leftarrow A_{2}D_{2}$ and $T\leftarrow
A_{3}D_{3}$.

19: return $C$, $R$ and $T$.

20:end procedure

######  Theorem C.21.

Given a 3rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any $k\geq
1$, there exists an algorithm which takes
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ time and outputs
three matrices: $C\in\mathbb{R}^{n\times c}$, a subset of columns of $A$,
$R\in\mathbb{R}^{n\times r}$ a subset of rows of $A$, and
$T\in\mathbb{R}^{n\times t}$, a subset of tubes of $A$, together with a tensor
$U\in\mathbb{R}^{c\times r\times t}$ with $\operatorname{rank}(U)=k^{\prime}$
where $c=r=t=\operatorname{poly}(k,1/\epsilon)$ and
$k^{\prime}=\operatorname{poly}(k,1/\epsilon)$ such that

$\displaystyle\|U(C,R,T)-A\|_{F}^{2}\leq(1+\epsilon)\underset{\operatorname{rank}-k\leavevmode\nobreak\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2},$

or equivalently,

$\displaystyle\left\|\sum_{i=1}^{c}\sum_{j=1}^{r}\sum_{l=1}^{t}U_{i,j,l}\cdot
C_{i}\otimes R_{j}\otimes
T_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\underset{\operatorname{rank}-k\leavevmode\nobreak\
A_{k}}{\min}\|A_{k}-A\|_{F}^{2}$

holds with probability $9/10$.

###### Proof.

The proof follows by combining Theorem 1.1 and Theorem 1.3 directly. ∎

###  C.7 CURT decomposition, $(1+\epsilon)$-approximation

####  C.7.1 Properties of leverage score sampling and BSS sampling

Notice that, the BSS algorithm is a deterministic procedure developed in
[BSS12] for selecting rows from a matrix $A\in\mathbb{R}^{n\times d}$ (with
$\|A\|_{2}\leq 1$ and $\|A\|_{F}^{2}\leq k$) using a selection matrix $S$ so
that

$\displaystyle\|A^{\top}S^{\top}SA-A^{\top}A\|_{2}\leq\epsilon.$

The algorithm runs in $\operatorname{poly}(n,d,1/\epsilon)$ time. Using the
ideas from [BW14] and [CEM+15], we are able to reduce the number of nonzero
entries from $O(\epsilon^{-2}k\log k)$ to $O(\epsilon^{-2}k)$, and also
improve the running time to input sparsity.

######  Lemma C.22 (Leverage score preserves subspace embedding - Theorem 2.11
in [Woo14]).

Given a $\operatorname{rank}$-$k$ matrix $A\in\mathbb{R}^{n\times d}$, via
leverage score sampling, we can obtain a diagonal matrix $D$ with $m$ nonzero
entries such that, letting $B=DA$, if $m=O(\epsilon^{-2}k\log k)$, then, with
probability at least $0.999$, for all $x\in\mathbb{R}^{d}$,

$\displaystyle(1-\epsilon)\|Ax\|_{2}\leq\|Bx\|_{2}\leq(1+\epsilon)\|Ax\|_{2}$

######  Lemma C.23.

Given a $\operatorname{rank}$-$k$ matrix $A\in\mathbb{R}^{n\times d}$, there
exists an algorithm that runs in
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ time and outputs
a matrix $B$ containing $O(\epsilon^{-2}k\log k)$ re-weighted rows of $A$,
such that with probability at least $0.999$, for all $x\in\mathbb{R}^{d}$,

$\displaystyle(1-\epsilon)\|Ax\|_{2}\leq\|Bx\|_{2}\leq(1+\epsilon)\|Ax\|_{2}$

###### Proof.

We choose a sparse embedding matrix (Definition B.16)
$\Pi\in\mathbb{R}^{d\times s}$ with $s=\operatorname{poly}(k/\epsilon)$. With
probability at least $0.999$, $\Pi^{\top}$ is a subspace embedding of
$A^{\top}$. Thus, $\operatorname{rank}(A\Pi)=\operatorname{rank}(A)$. Also,
the leverage scores of $A\Pi$ are the same as those of $A$. Thus, we can
compute the leverage scores of $A\Pi$. The running time of computing $A\Pi$ is
$O(\operatorname{nnz}(A))$. Thus the total running time is
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$. ∎

######  Lemma C.24.

Let $B$ denote a matrix which contains $O(\epsilon^{-2}k\log k)$ rows of
$A\in\mathbb{R}^{n\times d}$. Choosing $\Pi$ to be a sparse subspace embedding
matrix of size $d\times O(\epsilon^{-6}(k\log k)^{2})$, with probability at
least $0.999$,

$\displaystyle\|B\Pi\Pi^{\top}B^{\top}-BB^{\top}\|_{2}\leq\epsilon\|B\|_{2}^{2}.$

Combining Lemma C.23, C.24 and the BSS algorithm, we obtain:

######  Lemma C.25.

Given a $\operatorname{rank}$-$k$ matrix $A\in\mathbb{R}^{n\times d}$, there
exists an algorithm that runs in
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ time and outputs
a sampling and rescaling diagonal matrix $S$ that selects $O(\epsilon^{-2}k)$
re-weighted rows of $A$, such that, with probability at least $0.999$,

$\displaystyle\|A^{\top}S^{\top}SA-A^{\top}A\|_{2}\leq\epsilon\|A\|_{2}^{2}.$

or equivalently, for all $x\in\mathbb{R}^{d}$,

$\displaystyle(1-\epsilon)\|Ax\|_{2}\leq\|SAx\|_{2}\leq(1+\epsilon)\|Ax\|_{2}.$

###### Proof.

Using Lemma C.23, we can obtain $B$. Then we apply a sparse subspace embedding
matrix $\Pi$ on the right of $B$. At the end, we run the BSS algorithm on
$B\Pi$ and we are able to output $O(\epsilon^{-2}k)$ re-weighted rows of
$B\Pi$. Using these rows, we are able to determine $O(\epsilon^{-2}k)$ re-
weighted rows of $A$. ∎

####  C.7.2 Row sampling for linear regression

######  Theorem C.26 (Theorem 5 in [CNW15]).

We are given $A\in\mathbb{R}^{n\times d}$ with $\|A\|_{2}^{2}\leq 1$ and
$\|A\|_{F}^{2}\leq k$, and an $\epsilon\in(0,1)$. There exists a diagonal
matrix $S$ with $O(k/\epsilon^{2})$ nonzero entries such that

$\displaystyle\|(SA)^{\top}SA-A^{\top}A\|_{2}\leq\epsilon.$

######  Corollary C.27.

Given a rank-$k$ matrix $A\in\mathbb{R}^{n\times d}$, vector
$b\in\mathbb{R}^{n}$, and parameter $\epsilon>0$, let
$U\in\mathbb{R}^{n\times(k+1)}$ denote an orthonormal basis of $[A,b]$. Let
$S\in\mathbb{R}^{n\times n}$ denote a sampling and rescaling diagonal matrix
according to Leverage score sampling and sparse BSS sampling of $U$ with $m$
nonzero entries. If $m=O(k)$, then $S$ is a $(1\pm 1/2)$ subspace embedding
for $U$; if $m=O(k/\epsilon)$, then $S$ satisfies $\sqrt{\epsilon}$-operator
norm approximate matrix product for $U$.

###### Proof.

This follows by Lemma C.22, Lemma C.24 and Theorem C.26. ∎

######  Lemma C.28 ([NW14]).

Given $A\in\mathbb{R}^{n\times d}$ and $b\in\mathbb{R}^{n}$, let
$S\in\mathbb{R}^{n\times n}$ denote a sampling and rescaling diagonal matrix.
Let $x^{*}$ denote $\arg\min_{x}\|Ax-b\|_{2}^{2}$ and $x^{\prime}$ denote
$\arg\min_{x}\|SAx-Sb\|_{2}^{2}$. If $S$ is a $(1\pm 1/2)$ subspace embedding
for the column span of $A$, and $\epsilon^{\prime}$
(=$\sqrt{\epsilon}$)-operator norm approximate matrix product for $U$ adjoined
with $b-Ax^{*}$, then, with probability at least $.999$,

$\displaystyle\|Ax^{\prime}-b\|_{2}^{2}\leq(1+\epsilon)\|Ax^{*}-b\|_{2}^{2}.$

###### Proof.

We define $\operatorname{OPT}=\underset{x}{\min}\|Ax-b\|_{2}$. We define
$x^{\prime}=\underset{x}{\arg\min}\|SAx-Sb\|_{2}^{2}$ and
$x^{*}=\underset{x}{\arg\min}\|Ax-b\|_{2}^{2}$. Let $w=b-Ax^{*}$. Let $U$
denote an orthonormal basis of $A$. We can write $Ax^{\prime}-Ax^{*}=U\beta$.
Then, we have,

$\displaystyle\|Ax^{\prime}-b\|_{2}^{2}=$ $\displaystyle\leavevmode\nobreak\
\|Ax^{\prime}-Ax^{*}+AA^{\dagger}b-b\|_{2}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}x^{*}=A^{\dagger}b$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|U\beta+(UU^{\top}-I)b\|_{2}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|Ax^{*}-Ax^{\prime}\|_{2}^{2}+\|Ax^{*}-b\|_{2}^{2}$ by Pythagorean Theorem
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|U\beta\|_{2}^{2}+\operatorname{OPT}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|\beta\|_{2}^{2}+\operatorname{OPT}^{2}.$

If S is a $(1\pm 1/2)$ subspace embedding for $U$, then we can show

$\displaystyle\|\beta\|_{2}-\|U^{\top}S^{\top}SU\beta\|_{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\|\beta-U^{\top}S^{\top}SU\beta\|_{2}$ by triangle inequality $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|(I-U^{\top}S^{\top}SU)\beta\|_{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\|I-U^{\top}S^{\top}SU\|_{2}\cdot\|\beta\|_{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \frac{1}{2}\|\beta\|_{2}.$

Thus, we obtain

$\displaystyle\|U^{\top}S^{\top}SU\beta\|_{2}\geq\|\beta\|_{2}/2.$

Next, we can show

$\displaystyle\|U^{\top}S^{\top}SU\beta\|_{2}=$
$\displaystyle\leavevmode\nobreak\
\|U^{\top}S^{\top}S(Ax^{\prime}-Ax^{*})\|_{2}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|U^{\top}S^{\top}S(A(SA)^{\dagger}Sb-
Ax^{*})\|_{2}$ $\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}x^{\prime}=(SA)^{\dagger}Sb$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|U^{\top}S^{\top}S(b-Ax^{*})\|_{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}SA(SA)^{\dagger}=I$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|U^{\top}S^{\top}Sw\|_{2}.$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }w=b-Ax^{*}$

We define $U^{\prime}=\begin{bmatrix}U&w/\|w\|_{2}\end{bmatrix}$. We define
$X$ and $y$ to satisfy $U=U^{\prime}X$ and $w=U^{\prime}y$. Then, we have

$\displaystyle\leavevmode\nobreak\ \|U^{\top}S^{\top}Sw\|_{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|U^{\top}S^{\top}Sw-U^{\top}w\|_{2}$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }U^{\top}w=0$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|X^{\top}U^{\prime\top}S^{\top}SU^{\prime}y-X^{\top}U^{\prime\top}U^{\prime}y\|_{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|X^{\top}(U^{\prime\top}S^{\top}SU^{\prime}-I)y\|_{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\|X\|_{2}\cdot\|U^{\prime\top}S^{\top}SU^{\prime}-I\|_{2}\cdot\|y\|_{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\epsilon^{\prime}\|X\|_{2}\|y\|_{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \epsilon^{\prime}\|U\|_{2}\|w\|_{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\epsilon^{\prime}\operatorname{OPT},$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }\|U\|_{2}=1\text{\leavevmode\nobreak\
and\leavevmode\nobreak\ }\|w\|_{2}=\operatorname{OPT}$

where the fifth inequality follows since $S$ satisfies
$\epsilon^{\prime}$-operator norm approximate matrix product for the column
span of $U$ adjoined with $w$.

Putting it all together, we have

$\displaystyle\|Ax^{\prime}-b\|_{2}^{2}=$ $\displaystyle\leavevmode\nobreak\
\|Ax^{*}-b\|_{2}^{2}+\|Ax^{*}-Ax^{\prime}\|_{2}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \operatorname{OPT}^{2}+\|\beta\|_{2}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\operatorname{OPT}^{2}+4\|U^{\top}S^{\top}Sw\|_{2}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\operatorname{OPT}^{2}+4(\epsilon^{\prime}\operatorname{OPT})^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
(1+\epsilon)\operatorname{OPT}^{2}.$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }\epsilon^{\prime}=\frac{1}{2}\sqrt{\epsilon}.$

Finally, note that $S$ satisfies $\epsilon^{\prime}$-operator norm approximate
matrix product for $U$ adjoined with $w$ if it is a
$(1\pm\epsilon^{\prime})$-subspace embedding for $U$ adjoined with $w$, which
holds using BSS sampling by Theorem 5 of [CNW15] with $O(d/\epsilon)$ samples.
∎

####  C.7.3 Leverage scores for multiple regression

######  Lemma C.29 (see, e.g., Lemma 32 in [CW13] among other places).

Given matrix $A\in\mathbb{R}^{n\times d}$ with orthonormal columns, and
parameter $\epsilon>0$, if $S\in\mathbb{R}^{n\times n}$ is a sampling and
rescaling diagonal matrix according to the leverage scores of $A$ where the
number of nonzero entries is $t=O(1/\epsilon^{2})$, then, for any
$B\in\mathbb{R}^{n\times m}$, we have

$\displaystyle\|A^{\top}S^{\top}SB-A^{\top}B\|_{F}^{2}<\epsilon^{2}\|A\|_{F}^{2}\|B\|_{F}^{2},$

holds with probability at least $0.9999$.

######  Corollary C.30.

Given matrix $A\in\mathbb{R}^{n\times d}$ with orthonormal columns, and
parameter $\epsilon>0$, if $S\in\mathbb{R}^{n\times n}$ is a sampling and
rescaling diagonal matrix according to the leverage scores of $A$ with $m$
nonzero entries, then if $m=O(d\log d)$, then $S$ is a $(1\pm 1/2)$ subspace
embedding for $A$. If $m=O(d/\epsilon)$, then $S$ satisfies
$\sqrt{\epsilon/d}$-Frobenius norm approximate matrix product for $A$.

###### Proof.

This follows by Lemma C.22 and Lemma C.29. ∎

######  Lemma C.31 ([NW14]).

Given $A\in\mathbb{R}^{n\times d}$ and $B\in\mathbb{R}^{n\times m}$, let
$S\in\mathbb{R}^{n\times n}$ denote a sampling and rescaling matrix according
to $A$. Let $X^{*}$ denote $\arg\min_{X}\|AX-B\|_{F}^{2}$ and $X^{\prime}$
denote $\arg\min_{X}\|SAX-SB\|_{F}^{2}$. Let $U$ denote an orthonormal basis
for $A$. If $S$ is a $(1\pm 1/2)$ subspace embedding for $U$, and satisfies
$\epsilon^{\prime}$(=$\sqrt{\epsilon/d}$)-Frobenius norm approximate matrix
product for $U$, then, we have that

$\displaystyle\|AX^{\prime}-B\|_{F}^{2}\leq(1+\epsilon)\|AX^{*}-B\|_{F}^{2}$

holds with probability at least $0.999$.

###### Proof.

We define $\operatorname{OPT}=\min_{X}\|AX-B\|_{F}$. Let $A=U\Sigma V^{\top}$
denote the SVD of $A$. Since $A$ has rank $k$, $U$ and $V$ have $k$ columns.
We can write $A(X^{\prime}-X^{*})=U\beta$. Then, we have

$\displaystyle\|AX^{\prime}-B\|_{F}^{2}=$ $\displaystyle\leavevmode\nobreak\
\|AX^{\prime}-AX^{*}+AA^{\dagger}B-B\|_{F}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}X^{*}=A^{\dagger}B$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|U\beta+(UU^{\top}-I)B\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|AX^{*}-AX^{\prime}\|_{F}^{2}+\|AX^{*}-B\|_{F}^{2}$ by Pythagorean Theorem
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|U\beta\|_{F}^{2}+\operatorname{OPT}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|\beta\|_{F}^{2}+\operatorname{OPT}^{2}.$
(15)

If $S$ is a $(1\pm 1/2)$ subspace embedding for $U$, then we can show,

$\displaystyle\leavevmode\nobreak\
\|\beta\|_{F}-\|U^{\top}S^{\top}SSU\beta\|_{F}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \|\beta-U^{\top}S^{\top}SU\beta\|_{F}$ by
triangle inequality $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|(I-U^{\top}S^{\top}SU)\beta\|_{F}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\|(I-U^{\top}S^{\top}SU)\|_{2}\cdot\|\beta\|_{F}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}\|AB\|_{F}\leq\|A\|_{2}\|B\|_{F}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \frac{1}{2}\|\beta\|_{F}.$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}\|(I-U^{\top}S^{\top}SU)\|_{2}\leq 1/2$

Thus, we obtain

$\displaystyle\|U^{\top}S^{\top}SU\beta\|_{F}\geq\|\beta\|_{F}/2.$ (16)

Next, we can show

$\displaystyle\|U^{\top}S^{\top}SU\beta\|_{F}=$
$\displaystyle\leavevmode\nobreak\
\|U^{\top}S^{\top}S(AX^{\prime}-AX^{*})\|_{F}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|U^{\top}S^{\top}S(A(SA)^{\dagger}Sb-
AX^{*})\|_{F}$ $\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}X^{\prime}=(SA)^{\dagger}SB$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \|U^{\top}S^{\top}S(B-AX^{*})\|_{F}.$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}SA(SA)^{\dagger}=I$

Then, we can show

$\displaystyle\|U^{\top}S^{\top}S(B-AX^{*})\|_{F}\leq$
$\displaystyle\leavevmode\nobreak\
\epsilon^{\prime}\|U^{\top}\|_{F}\|B-AX^{*}\|_{F}$ by Lemma C.29
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\epsilon^{\prime}\sqrt{d}\operatorname{OPT}.$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}\|U\|_{F}=\sqrt{d}\text{\leavevmode\nobreak\ and\leavevmode\nobreak\
}\|B-AX^{*}\|_{F}=\operatorname{OPT}$ (17)

Putting it all together, we have

$\displaystyle\|AX^{\prime}-B\|_{F}^{2}=$ $\displaystyle\leavevmode\nobreak\
\|AX^{*}-B\|_{F}^{2}+\|AX^{*}-AX^{\prime}\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \operatorname{OPT}^{2}+\|\beta\|_{F}^{2}$
by Equation (C.7.3) $\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\operatorname{OPT}^{2}+4\|U^{\top}S^{\top}Sw\|_{F}^{2}$ by Equation (16)
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\operatorname{OPT}^{2}+4(\epsilon^{\prime}\sqrt{d}\operatorname{OPT})^{2}$ by
Equation (C.7.3) $\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
(1+\epsilon)\operatorname{OPT}^{2}.$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }\epsilon^{\prime}=\frac{1}{2}\sqrt{\epsilon/d}$

∎

####  C.7.4 Sampling columns according to leverage scores implicitly,
improving polynomial running time to nearly linear running time

This section explains an algorithm that is able to sample from the leverage
scores from the $\odot$ product of two matrices $U,V$ without explicitly
writing down $U\odot V$. To build this algorithm we combine TensorSketch, some
ideas from [DMIMW12] and some ideas from [AKO11, MW10]. Finally, we are able
to improve the running time of sampling columns according to leverage scores
from $\Omega(n^{2})$ to $\widetilde{O}(n)$. Given two matrices
$U,V\in\mathbb{R}^{k\times n}$, we define $A\in\mathbb{R}^{k\times
n_{1}n_{2}}$ to be the matrix where the $i$-th row of $A$ is the vectorization
of $U^{i}\otimes V^{i}$, $\forall i\in[k]$. Naïvely, in order to sample
$O(\operatorname{poly}(k,1/\epsilon))$ rows from $A^{\top}$ according to
leverage scores, we need to write down $n^{2}$ leverage scores. This approach
will take at least $\Omega(n^{2})$ running time. In the rest of this section,
we will explain how to do it in $O(n\cdot\operatorname{poly}(\log
n,k,1/\epsilon))$ time. In Section C.10.1, we will explain how to extend this
idea from $3$rd order tensors to general $q$-th order tensors and remove the
$\operatorname{poly}(\log n)$ factor from running time, i.e., obtain
$O(n\cdot\operatorname{poly}(k,1/\epsilon))$ time.

Algorithm 10 Fast Tensor Leverage Score Sampling

1:procedure
FastTensorLeverageScore($U,V,n_{1},n_{2},k,\epsilon,R_{\text{samples}}$)
$\triangleright$ Lemma C.32

2: $s_{1}\leftarrow\operatorname{poly}(k,1/\epsilon)$.

3: $g_{1}\leftarrow g_{2}\leftarrow g_{3}\leftarrow
O(\epsilon^{-2}\log(n_{1}n_{2}))$.

4: Choose $\Pi\in\mathbb{R}^{n_{1}n_{2}\times s_{1}}$ to be a TensorSketch.
$\triangleright$ Definition B.34

5: Compute $R^{-1}\in\mathbb{R}^{k\times k}$ by using $(U\odot V)\Pi$.
$\triangleright$ $U\in\mathbb{R}^{k\times n_{1}}$, $V\in\mathbb{R}^{k\times
n_{2}}$

6: Choose $G_{1}\in\mathbb{R}^{g_{1}\times k}$ to be a Gaussian sketching
matrix.

7: for $i=1\to g_{1}$ do

8: $w\leftarrow(G^{i}R^{-1})^{\top}$ $\triangleright$ $G^{i}$ denotes the
$i$-th row of $G$

9: for $j=1\to[n_{1}]$ do $\triangleright$ Form matrix $U^{\prime
i}\in\mathbb{R}^{k\times n_{1}}$

10: $U^{\prime i}_{j}\leftarrow w\circ U_{j},\forall j\in[n_{1}]$.
$\triangleright$ $U_{j}$ denotes the $j$-th column of $U\in\mathbb{R}^{k\times
n_{1}}$

11: end for

12: end for

13: Choose $G_{2,i}\in\mathbb{R}^{g_{2}\times n_{1}}$ to be a Gaussian
sketching matrix.

14: for $i=1\to g_{1}$ do

15: $\alpha_{i}\leftarrow\|(G_{2,i}U^{\prime i\top})V\|_{F}^{2}$.

16: Choose $G_{3,i}\in\mathbb{R}^{g_{3}\times n_{1}}$ to be a Gaussian
sketching matrix.

17: for $j_{2}=1\to n_{2}$ do

18: $\beta_{i,j}\leftarrow\|G_{3,i}(U^{\prime i\top})V_{j_{2}}\|_{2}^{2}$.

19: end for

20: end for

21: ${\cal S}\leftarrow\emptyset$.

22: for $r=1\to R_{\text{samples}}$ do

23: Sample $i$ from $[g_{1}]$ with probability
$\alpha_{i}/\sum_{i^{\prime}=1}^{g_{1}}\alpha_{i^{\prime}}$.

24: Sample $j_{2}$ from $[n_{2}]$ with probability
$\beta_{i,j_{2}}/\sum_{j^{\prime}_{2}=1}^{n_{2}}\beta_{i,j^{\prime}_{2}}$.

25: for $j_{1}=1\to n_{1}$ do

26: $\gamma_{j_{1}}\leftarrow((U^{\prime i\top})^{j_{1}}V_{j_{2}})^{2}$.

27: end for

28: Sample $j_{1}$ from $[n_{1}]$ with probability
$\gamma_{j_{1}}/\sum_{j^{\prime}_{1}=1}^{n_{1}}\gamma_{j^{\prime}_{1}}$.

29: ${\cal S}\leftarrow{\cal S}\cup(j_{1},j_{2})$.

30: end for

31: Convert ${\cal S}$ into a diagonal matrix $D$ with at most
$R_{\text{samples}}$ nonzero entries.

32: return $D$. $\triangleright$ Diagonal matrix
$D\in\mathbb{R}^{n_{1}n_{2}\times n_{1}n_{2}}$

33:end procedure

######  Lemma C.32.

Given two matrices $U\in\mathbb{R}^{k\times n_{1}}$ and
$V\in\mathbb{R}^{k\times n_{2}}$, there exists an algorithm that takes
$O((n_{1}+n_{2})\cdot\operatorname{poly}(\log(n_{1}n_{2}),k)\cdot
R_{\mathrm{samples}})$ time and samples $R_{\mathrm{samples}}$ columns of
$U\odot V\in\mathbb{R}^{k\times n_{1}n_{2}}$ according to the leverage scores
of $R^{-1}(U\odot V)$, where $R$ is the $R$ of a QR factorization.

###### Proof.

We choose $\Pi\in\mathbb{R}^{n_{1}n_{2}\times s_{1}}$ to be a TensorSketch.
Then, according to Section B.10, we can compute $R^{-1}$ in
$n\cdot\operatorname{poly}(\log n,k,1/\epsilon)$ time, where $R$ is the $R$ in
a QR-factorization. We want to sample columns from $U\odot V$ according to the
square of the $\ell_{2}$-norms of each column of $R^{-1}(U\odot V)$. However,
explicitly writing down the matrix $R^{-1}(U\odot V)$ takes $kn_{1}n_{2}$
time, and the number of columns is already $n_{1}n_{2}$. The goal is to sample
columns from $R^{-1}(U\odot V)$ without explicitly computing the square of the
$\ell_{2}$-norm of each column.

The first simple observation is that the following two sampling procedures are
equivalent in terms of the column samples of a matrix that they take. (1) We
sample a single entry from the matrix $R^{-1}(U\odot V)$ proportional to its
squared value. (2) We sample a column from the matrix $R^{-1}(U\odot V)$
proportional to its squared $\ell_{2}$-norm. Let the $(i,j_{1},j_{2})$-th
entry denote the entry in the $i$-th row and the $(j_{1}-1)n_{2}+j_{2}$-th
column. We can show, for a particular column $(j_{1}-1)n_{2}+j_{2}$,

$\displaystyle\leavevmode\nobreak\ \Pr[\text{sample\leavevmode\nobreak\
an\leavevmode\nobreak\ entry\leavevmode\nobreak\ from\leavevmode\nobreak\
the\leavevmode\nobreak\ }(j_{1}-1)n_{2}+j_{2}\text{\leavevmode\nobreak\
th\leavevmode\nobreak\ column\leavevmode\nobreak\ of\leavevmode\nobreak\
a\leavevmode\nobreak\ matrix}]$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\sum_{i=1}^{k}\Pr[\text{sample\leavevmode\nobreak\ the}\leavevmode\nobreak\
(i,j_{1},j_{2})\text{-th\leavevmode\nobreak\ entry\leavevmode\nobreak\
of\leavevmode\nobreak\ matrix}]$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \sum_{i=1}^{k}\frac{|(R^{-1}(U\odot
V))_{i,(j_{1}-1)n_{2}+j_{2}}|^{2}}{\|R^{-1}(U\odot V)\|_{F}^{2}}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\ \frac{\|(R^{-1}(U\odot
V))_{(j_{1}-1)n_{2}+j_{2}}\|^{2}}{\|R^{-1}(U\odot V)\|_{F}^{2}}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\Pr[\text{sample\leavevmode\nobreak\ the\leavevmode\nobreak\
}(j_{1}-1)n_{2}+j_{2}\text{\leavevmode\nobreak\ th\leavevmode\nobreak\
column\leavevmode\nobreak\ of\leavevmode\nobreak\ matrix}].$ (18)

Thus, it is sufficient to show how to sample a single entry from matrix
$R^{-1}(U\odot V)$ proportional to its squared value without writing down all
of the entries of a $k\times n_{1}n_{2}$ matrix.

We choose a Gaussian matrix $G_{1}\in\mathbb{R}^{g_{1}\times k}$ with
$g_{1}=O(\epsilon^{-2}\log(n_{1}n_{2}))$. By Claim C.33 we can reduce the
length of each column vector of matrix $R^{-1}U\odot V$ from $k$ to $g_{1}$
while preserving the squared $\ell_{2}$-norm of all columns simultaneously.
Thus, we obtain a new matrix $GR^{-1}(U\odot V)\in\mathbb{R}^{g_{1}\times
n_{1}n_{2}}$, and sampling from this new matrix is equivalent to sampling from
the original matrix $R^{-1}(U\odot V)$.

In the following paragraphs, we explain a sampling procedure (also described
in Procedure FastTensorLeverageScore in Algorithm 10) which contains three
sampling steps. The first step is sampling $i$ from $[g_{1}]$, the second step
is sampling $j_{2}$ from $[n_{2}]$, and the last step is sampling $j_{1}$ from
$[n_{1}]$.

For each $j_{1}\in[n_{1}]$, let $U_{j_{1}}$ denote the $j_{1}$-th column of
$U$. For each $i\in[g_{1}]$, let $G_{1}^{i}$ denote the $i$-th row of matrix
$G_{1}\in\mathbb{R}^{g_{1}\times k}$, let $U^{\prime i}\in\mathbb{R}^{k\times
n_{1}}$ denote a matrix where the $j_{1}$-th column is
$(G^{i}R^{-1})^{\top}\circ U_{j_{1}}\in\mathbb{R}^{k},\forall j\in[n_{1}]$.
Then, using Claim C.37, we have that $(G^{i}R^{-1})\cdot(U\odot
V)\in\mathbb{R}^{n_{1}n_{2}}$ is a row vector where the entry in the
$(j_{1}-1)n_{2}+j_{2}$-th coordinate is the entry in the $j_{1}$-th row and
$j_{2}$-th column of matrix $(U^{\prime i\top}V)\in\mathbb{R}^{n_{1}\times
n_{2}}$. Further, the squared $\ell_{2}$-norm of vector
$(G^{i}R^{-1})\cdot(U\odot V)$ is equal to the squared Frobenius norm of
matrix $(U^{\prime i\top}V)$. Thus, sampling $i$ proportional to the squared
$\ell_{2}$-norm of vector $(G^{i}R^{-1})\cdot(U\odot V)$ is equivalent to
sampling $i$ proportional to the squared Frobenius norm of matrix $(U^{\prime
i\top}V)$. Naïvely, computing the Frobenius norm of an $n_{1}\times n_{2}$
matrix requires $O(n_{1}n_{2})$ time. However, we can choose a Gaussian matrix
$G_{2,i}\in\mathbb{R}^{g_{2}\times n_{1}}$ to sample according to the value
$\|(G_{2,i}U^{\prime i\top})V\|_{F}^{2}$, which can be computed in
$O((n_{1}+n_{2})g_{2}k)$ time. By claim C.35, $\|(G_{2,i}U^{\prime
i\top})V\|_{F}^{2}\approx\|(U^{\prime i\top})V\|_{F}^{2}$ with high
probability. So far, we have finished the first step of the sampling
procedure.

For the second step of the sampling procedure, we need to sample $j_{2}$ from
$[n_{2}]$. To do that, we need to compute the squared $\ell_{2}$-norm of each
column of $U^{\prime i\top}V\in\mathbb{R}^{n_{1}\times n_{2}}$. This can be
done by choosing another Gaussian matrix $G_{3,i}\in\mathbb{R}^{g_{3}\times
n_{1}}$. For all $j_{2}\in[n_{2}]$, by Claim C.36, we have $\|G_{3,i}U^{\prime
i\top}V_{j_{2}}\|_{2}^{2}\approx\|U^{\prime i\top}V_{j_{2}}\|_{2}^{2}$. Also,
for $j_{2}\in[n_{2}]$, $\|G_{3,i}U^{\prime i\top}V_{j_{2}}\|_{2}^{2}$ can be
computed in nearly linear in $n_{1}+n_{2}$ time.

For the third step of the sampling procedure, we need to sample $j_{1}$ from
$[n_{1}]$. Since we already have $i$ and $j_{2}$ from the previous two steps,
we can directly compute $|(U^{\prime i\top})^{j_{1}}V_{j_{2}}|^{2}$, for all
$j_{1}$. This only takes $O(n_{1}k)$ time.

Overall, the running time is
$O((n_{1}+n_{2})\cdot\operatorname{poly}(\log(n_{1}n_{2}),k,1/\epsilon))$.
Because our estimates are accurate enough, our sampling probabilities are also
good approximations to the leverage score sampling probabilities. Putting it
all together, we complete the proof. ∎

######  Claim C.33.

Given matrix $R^{-1}(U\odot V)\in\mathbb{R}^{k\times n_{1}n_{2}}$, let
$G_{1}\in\mathbb{R}^{g_{1}\times k}$ denote a Gaussian matrix with
$g_{1}=(\epsilon^{-2}\log(n_{1}n_{2}))$. Then with probability at least
$1-1/\operatorname{poly}(n_{1}n_{2})$, we have: for all $j\in[n_{1}n_{2}]$,

$\displaystyle(1-\epsilon)\|R^{-1}(U\odot
V)_{j}\|_{2}^{2}\leq\|G_{1}R^{-1}(U\odot
V)_{j}\|_{2}^{2}\leq(1+\epsilon)\|R^{-1}(U\odot V)_{j}\|_{2}^{2}.$

###### Proof.

This follows by the Johnson-Lindenstrauss Lemma. ∎

######  Claim C.34.

For a fixed $i\in[g_{1}]$, let $G_{2,i}\in\mathbb{R}^{g_{2}\times n_{1}}$
denote a Gaussian matrix with $g_{2}=O(\epsilon^{-2}\log(n_{1}n_{2}))$. Then
with probability at least $1-1/\operatorname{poly}(n_{1}n_{2})$, we have: for
all $j_{2}\in[n_{2}]$,

$\displaystyle(1-\epsilon)\|U^{\prime
i\top}V_{j_{2}}\|_{2}^{2}\leq\|(G_{2,i}U^{\prime
i\top})V_{j_{2}}\|_{2}\leq(1+\epsilon)\|U^{\prime i\top}V_{j_{2}}\|_{2}^{2}.$

By taking the union bound over all $i\in[g_{1}]$, we obtain a stronger claim,

######  Claim C.35.

With probability at least $1-1/\operatorname{poly}(n_{1}n_{2})$, we have : for
all $i\in[g_{1}]$, for all $j_{2}\in[n_{2}]$,

$\displaystyle(1-\epsilon)\|U^{\prime
i\top}V_{j_{2}}\|_{2}^{2}\leq\|(G_{2,i}U^{\prime
i\top})V_{j_{2}}\|_{2}\leq(1+\epsilon)\|U^{\prime i\top}V_{j_{2}}\|_{2}^{2}.$

Similarly, if we choose $G_{3,i}$ to be a Gaussian matrix, we can obtain the
same result as for $G_{2,i}$:

######  Claim C.36.

With probability at least $1-1/\operatorname{poly}(n_{1}n_{2})$, we have : for
all $i\in[g_{1}]$, for all $j_{2}\in[n_{2}]$,

$\displaystyle(1-\epsilon)\|U^{\prime
i\top}V_{j_{2}}\|_{2}^{2}\leq\|(G_{3,i}U^{\prime
i\top})V_{j_{2}}\|_{2}\leq(1+\epsilon)\|U^{\prime i\top}V_{j_{2}}\|_{2}^{2}.$

######  Claim C.37.

For any $i\in[g_{1}]$, $j_{1}\in[n_{1}]$, $j_{2}\in[n_{2}]$, let $G_{1}^{i}$
denote the $i$-th row of matrix $G_{1}\in\mathbb{R}^{g_{1}\times k}$. Let
$(U\odot V)_{(j_{1}-1)n_{2}+j_{2}}$ denote the $(j_{1}-1)n_{2}+j_{2}$-th
column of matrix $\mathbb{R}^{k\times n_{1}n_{2}}$. Let $(U^{\prime
i\top})^{j_{1}}$ denote the $j_{1}$-th row of matrix $(U^{\prime
i\top})\in\mathbb{R}^{n_{1}\times k}$. Let $V_{j_{2}}$ denote the $j_{2}$-th
column of matrix $V\in\mathbb{R}^{k\times n_{2}}$. Then, we have

$\displaystyle G_{1}^{i}R^{-1}(U\odot V)_{(j_{1}-1)n_{2}+j_{2}}=(U^{\prime
i\top})^{j_{1}}V_{j_{2}}.$

###### Proof.

This follows by,

$\displaystyle G_{1}^{i}R^{-1}(U\odot
V)_{(j_{1}-1)n_{2}+j_{2}}=G_{1}^{i}R^{-1}(U_{j_{1}}\circ
V_{j_{2}})=(G_{1}^{i}R^{-1}\circ(U_{j_{1}})^{\top})V_{j_{2}}=(U^{\prime
i\top})^{j_{1}}V_{j_{2}}.$

∎

######  Lemma C.38.

Given $A\in\mathbb{R}^{n\times n^{2}}$, $V,W\in\mathbb{R}^{k\times n}$, for
any $\epsilon>0$, there exists an algorithm that runs in
$O(n\cdot\operatorname{poly}(k,1/\epsilon))$ time and outputs a diagonal
matrix $D\in\mathbb{R}^{n^{2}\times n^{2}}$ with $m=O(k\log k+k/\epsilon)$
nonzero entries such that,

$\displaystyle\|\widehat{U}(V\odot
W)-A\|_{F}^{2}\leq(1+\epsilon)\min_{U\in\mathbb{R}^{n\times k}}\|U(V\odot
W)-A\|_{F}^{2},$

holds with probability at least $0.999$, where $\widehat{U}$ denotes the
optimal solution to $\min_{U}\|U(V\odot W)D-AD\|_{F}^{2}$.

###### Proof.

This follows by combining Theorem C.46, Corollary C.30, and Lemma C.31. ∎

######  Remark C.39.

Replacing Theorem C.46 (Algorithm 15) by Lemma C.32 (Algorithm 10), we can
obtain a slightly different version of Lemma C.38 with
$n\operatorname{poly}(\log n,k,1/\epsilon)$ running time, where the dependence
on $k$ is better.

####  C.7.5 Input sparsity time algorithm

Algorithm 11 Frobenius Norm CURT Decomposition Algorithm, Input Sparsity Time
and Nearly Optimal Number of Samples

1:procedure FCURTInputSparsity($A,U_{B},V_{B},W_{B},n,k,\epsilon$)
$\triangleright$ Theorem C.40

2: $d_{1}\leftarrow d_{2}\leftarrow d_{3}\leftarrow O(k\log k+k/\epsilon)$.

3: $\epsilon_{0}\leftarrow 0.01$.

4: Form $B_{1}=V_{B}^{\top}\odot W_{B}^{\top}\in\mathbb{R}^{k\times n^{2}}$.

5:
$D_{1}\leftarrow$FastTensorLeverageScoreGeneralOrder($V_{B}^{\top},W_{B}^{\top},n,n,k,\epsilon_{0},d_{1}$).
$\triangleright$ Algorithm 15

6: Form $\widehat{U}=A_{1}D_{1}(B_{1}D_{1})^{\dagger}\in\mathbb{R}^{n\times
k}$.

7: Form $B_{2}=\widehat{U}^{\top}\odot W_{B}^{\top}\in\mathbb{R}^{k\times
n^{2}}$.

8:
$D_{2}\leftarrow$FastTensorLeverageScoreGeneralOrder($\widehat{U}^{\top},W_{B}^{\top},n,n,k,\epsilon_{0},d_{2}$).

9: Form $\widehat{V}=A_{2}D_{2}(B_{2}D_{2})^{\dagger}\in\mathbb{R}^{n\times
k}$.

10: Form
$B_{3}=\widehat{U}^{\top}\odot\widehat{V}^{\top}\in\mathbb{R}^{k\times
n^{2}}$.

11:
$D_{3}\leftarrow$FastTensorLeverageScoreGeneralOrder($\widehat{U}^{\top},\widehat{V}^{\top},n,n,k,\epsilon_{0},d_{3}$).

12: $C\leftarrow A_{1}D_{1}$, $R\leftarrow A_{2}D_{2}$, $T\leftarrow
A_{3}D_{3}$.

13:
$U\leftarrow\sum_{i=1}^{k}((B_{1}D_{1})^{\dagger})_{i}\otimes((B_{2}D_{2})^{\dagger})_{i}\otimes((B_{3}D_{3})^{\dagger})_{i}$.

14: return $C$, $R$, $T$ and $U$.

15:end procedure

######  Theorem C.40.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, let $k\geq
1$, and let $U_{B},V_{B},W_{B}\in\mathbb{R}^{n\times k}$ denote a rank-$k$,
$\alpha$-approximation to $A$. Then there exists an algorithm which takes
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ time and outputs
three matrices $C\in\mathbb{R}^{n\times c}$ with columns from $A$,
$R\in\mathbb{R}^{n\times r}$ with rows from $A$, $T\in\mathbb{R}^{n\times t}$
with tubes from $A$, and a tensor $U\in\mathbb{R}^{c\times r\times t}$ with
$\operatorname{rank}(U)=k$ such that $c=r=t=O(k\log k+k/\epsilon)$, and

$\displaystyle\left\|\sum_{i=1}^{c}\sum_{j=1}^{r}\sum_{l=1}^{t}U_{i,j,l}\cdot
C_{i}\otimes R_{j}\otimes
T_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\alpha\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}$

holds with probability $9/10$.

###### Proof.

We define

$\displaystyle\operatorname{OPT}:=\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}.$

We already have three matrices $U_{B}\in\mathbb{R}^{n\times k}$,
$V_{B}\in\mathbb{R}^{n\times k}$ and $W_{B}\in\mathbb{R}^{n\times k}$ and
these three matrices provide a $\operatorname{rank}$-$k$,
$\alpha$-approximation to $A$, i.e.,

$\displaystyle\left\|\sum_{i=1}^{k}(U_{B})_{i}\otimes(V_{B})_{i}\otimes(W_{B})_{i}-A\right\|_{F}^{2}\leq\alpha\operatorname{OPT}.$
(19)

Let $B_{1}=V_{B}^{\top}\odot W_{B}^{\top}\in\mathbb{R}^{k\times n^{2}}$ denote
the matrix where the $i$-th row is the vectorization of
$(V_{B})_{i}\otimes(W_{B})_{i}$. Let $D_{1}\in\mathbb{R}^{n^{2}\times n^{2}}$
be a sampling and rescaling matrix corresponding to sampling by the leverage
scores of $B_{1}^{\top}$; there are $d_{1}$ nonzero entries on the diagonal of
$D_{1}$. Let $A_{i}\in\mathbb{R}^{n\times n^{2}}$ denote the matrix obtained
by flattening $A$ along the $i$-th direction, for each $i\in[3]$.

Define $U^{*}\in\mathbb{R}^{n\times k}$ to be the optimal solution to
$\underset{U\in\mathbb{R}^{n\times k}}{\min}\|UB_{1}-A_{1}\|_{F}^{2}$,
$\widehat{U}=A_{1}D_{1}(B_{1}D_{1})^{\dagger}\in\mathbb{R}^{n\times k}$, and
$V_{0}\in\mathbb{R}^{n\times k}$ to be the optimal solution to
$\underset{V\in\mathbb{R}^{n\times k}}{\min}\|V\cdot(\widehat{U}^{\top}\odot
W_{B}^{\top})-A_{2}\|_{F}^{2}$. Due to Lemma C.38, if $d_{1}=O(k\log
k+k/\epsilon)$ then with constant probability, we have

$\displaystyle\|\widehat{U}B_{1}-A_{1}\|_{F}^{2}\leq\alpha_{D_{1}}\|U^{*}B_{1}-A_{1}\|_{F}^{2}.$
(20)

Recall that $(\widehat{U}^{\top}\odot W_{B}^{\top})\in\mathbb{R}^{k\times
n^{2}}$ denotes the matrix where the $i$-th row is the vectorization of
$\widehat{U}_{i}\otimes(W_{B})_{i}$, $\forall i\in[k]$. Now, we can show,

$\displaystyle\|V_{0}\cdot(\widehat{U}^{\top}\odot
W_{B}^{\top})-A_{2}\|_{F}^{2}\leq$ $\displaystyle\leavevmode\nobreak\
\|\widehat{U}B_{1}-A_{1}\|_{F}^{2}$ $\displaystyle\text{\leavevmode\nobreak\
by\leavevmode\nobreak\ }V_{0}=\underset{V\in\mathbb{R}^{n\times
k}}{\arg\min}\|V\cdot(\widehat{U}^{\top}\odot W_{B}^{\top})-A_{2}\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\alpha_{D_{1}}\|U^{*}B_{1}-A_{1}\|_{F}^{2}$ by Equation (20)
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\alpha_{D_{1}}\|U_{B}B_{1}-A_{1}\|_{F}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}U^{*}=\underset{U\in\mathbb{R}^{n\times
k}}{\arg\min}\|UB_{1}-A_{1}\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \alpha_{D_{1}}\alpha\operatorname{OPT}.$ by
Equation (19) (21)

We define $B_{2}=\widehat{U}^{\top}\odot W_{B}^{\top}$. Let
$D_{2}\in\mathbb{R}^{n^{2}\times n^{2}}$ be a sampling and rescaling matrix
corresponding to the leverage scores of $B_{2}^{\top}$. Suppose there are
$d_{2}$ nonzero entries on the diagonal of $D_{2}$.

Define $V^{*}\in\mathbb{R}^{n\times k}$ to be the optimal solution to
$\min_{V\in\mathbb{R}^{n\times k}}\|VB_{2}-A_{2}\|_{F}^{2}$,
$\widehat{V}=A_{2}D_{2}(B_{2}D_{2})^{\dagger}\in\mathbb{R}^{n\times k}$,
$W_{0}\in\mathbb{R}^{n\times k}$ to be the optimal solution to
$\underset{W\in\mathbb{R}^{n\times
k}}{\min}\|W\cdot(\widehat{U}^{\top}\odot\widehat{V}^{\top})-A_{3}\|_{F}^{2}$,
and $V^{\prime}$ to be the optimal solution to
$\underset{V\in\mathbb{R}^{n\times
k}}{\min}\|VB_{2}D_{2}-A_{2}D_{2}\|_{F}^{2}$.

Due to Lemma C.38, with constant probability, we have

$\displaystyle\|\widehat{V}B_{2}-A_{2}\|_{F}^{2}\leq\alpha_{D_{2}}\|V^{*}B_{2}-A_{2}\|_{F}^{2}.$
(22)

Recall that $(\widehat{U}^{\top}\odot\widehat{V}^{\top})\in\mathbb{R}^{k\times
n^{2}}$ denotes the matrix where the $i$-th row is the vectorization of
$\widehat{U}_{i}\otimes\widehat{V}_{i}$, $\forall i\in[k]$. Now, we can show,

$\displaystyle\|W_{0}\cdot(\widehat{U}^{\top}\odot\widehat{V}^{\top})-A_{3}\|_{F}^{2}\leq$
$\displaystyle\leavevmode\nobreak\ \|\widehat{V}B_{2}-A_{2}\|_{F}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}W_{0}=\underset{W\in\mathbb{R}^{n\times
k}}{\arg\min}\|W\cdot(\widehat{U}^{\top}\odot\widehat{V}^{\top})-A_{3}\|_{F}^{2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\alpha_{D_{2}}\|V^{*}B_{2}-A_{2}\|_{F}^{2}$ by Equation (22)
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\alpha_{D_{2}}\|V_{0}B_{2}-A_{2}\|_{F}^{2}$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}V^{*}=\underset{V\in\mathbb{R}^{n\times
k}}{\arg\min}\|VB_{2}-A_{2}\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\alpha_{D_{2}}\alpha_{D_{1}}\alpha\operatorname{OPT}.$ by Equation (C.7.5)
(23)

We define $B_{3}=\widehat{U}^{\top}\odot\widehat{V}^{\top}$. Let
$D_{3}\in\mathbb{R}^{n^{2}\times n^{2}}$ denote a sampling and rescaling
matrix corresponding to sampling by the leverage scores of $B_{3}^{\top}$.
Suppose there are $d_{3}$ nonzero entries on the diagonal of $D_{3}$.

Define $W^{*}\in\mathbb{R}^{n\times k}$ to be the optimal solution to
$\min_{W\in\mathbb{R}^{n\times k}}\|WB_{3}-A_{3}\|_{F}^{2}$,
$\widehat{W}=A_{3}D_{3}(B_{3}D_{3})^{\dagger}\in\mathbb{R}^{n\times k}$, and
$W^{\prime}$ to be the optimal solution to $\underset{W\in\mathbb{R}^{n\times
k}}{\min}\|WB_{3}D_{3}-A_{3}D_{3}\|_{F}^{2}$.

Due to Lemma C.38 with constant probability, we have

$\displaystyle\|\widehat{W}B_{3}-A_{3}\|_{F}^{2}\leq\alpha_{D_{3}}\|W^{*}B_{3}-A_{3}\|_{F}^{2}.$
(24)

Now we can show,

$\displaystyle\|\widehat{W}B_{3}-A_{3}\|_{F}^{2}\leq$
$\displaystyle\leavevmode\nobreak\
\alpha_{D_{3}}\|W^{*}B_{3}-A_{3}\|_{F}^{2},$ by Equation (24)
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\alpha_{D_{3}}\|W_{0}B_{3}-A_{3}\|_{F}^{2},$
$\displaystyle\text{\leavevmode\nobreak\ by\leavevmode\nobreak\
}W^{*}=\underset{W\in\mathbb{R}^{n\times
k}}{\arg\min}\|WB_{3}-A_{3}\|_{F}^{2}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\alpha_{D_{3}}\alpha_{D_{2}}\alpha_{D_{1}}\alpha\operatorname{OPT}.$ by
Equation (C.7.5)

This implies,

$\displaystyle\left\|\sum_{i=1}^{k}\widehat{U}_{i}\otimes\widehat{V}_{i}\otimes\widehat{W}_{i}-A\right\|_{F}^{2}\leq
O(1)\alpha\operatorname{OPT}^{2}.$

where $\widehat{U}=A_{1}D_{1}(B_{1}D_{1})^{\dagger}$,
$\widehat{V}=A_{2}D_{2}(B_{2}D_{2})^{\dagger}$,
$\widehat{W}=A_{3}D_{3}(B_{3}D_{3})^{\dagger}$.

By Lemma C.38, we need to set $d_{1}=d_{2}=d_{3}=O(k\log k+k/\epsilon)$. Note
that $B_{1}=(V_{B}^{\top}\odot W_{B}^{\top})$. Thus $D_{1}$ can be found in
$n\cdot\operatorname{poly}(k,1/\epsilon)$ time. Because $D_{1}$ has a small
number of nonzero entries on the diagonal, we can compute $B_{1}D_{1}$ quickly
without explicitly writing down $B_{1}$. Also $A_{1}D_{1}$ can be computed in
$\operatorname{nnz}(A)$ time. Using $(A_{1}D_{1})$ and $(B_{1}D_{1})$, we can
compute $\widehat{U}$ in $n\operatorname{poly}(k,1/\epsilon)$ time. In a
similar way, we can compute $B_{2}$, $D_{2}$, $B_{3}$, and $D_{3}$. Since
tensor $U$ is constructed based on three $\operatorname{poly}(k,1/\epsilon)$
size matrices, $(B_{1}D_{1})^{\dagger}$, $(B_{2}D_{2})^{\dagger}$, and
$(B_{3}D_{3})^{\dagger}$, the overall running time is
$O(\operatorname{nnz}(A)+n\operatorname{poly}(k,1/\epsilon))$ ∎

####  C.7.6 Optimal sample complexity algorithm

Algorithm 12 Frobenius Norm CURT Decomposition Algorithm, Optimal Sample
Complexity

1:procedure FCURTOptimalSamples($A,U_{B},V_{B},W_{B},n,k$) $\triangleright$
Theorem C.41

2: $d_{1}\leftarrow d_{2}\leftarrow d_{3}\leftarrow O(k/\epsilon)$.

3: Form $B_{1}=V_{B}^{\top}\odot W_{B}^{\top}\in\mathbb{R}^{k\times n^{2}}$.

4:
$D_{1}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{1}^{\top},B_{1}^{\top},n^{2},n,k,\epsilon$).
$\triangleright$ Algorithm 7

5: Let $d_{1}$ denote the number of nonzero entries in $D_{1}$.
$\triangleright$ $d_{1}=O(k/\epsilon)$

6: Form $\widehat{U}=A_{1}D_{1}(B_{1}D_{1})^{\dagger}\in\mathbb{R}^{n\times
k}$.

7: Form $B_{2}=\widehat{U}^{\top}\odot W_{B}^{\top}\in\mathbb{R}^{k\times
n^{2}}$.

8:
$D_{2}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{2}^{\top},B_{2}^{\top},n^{2},n,k,\epsilon$).
$\triangleright$ Algorithm 7

9: Let $d_{2}$ denote the number of nonzero entries in $D_{2}$.
$\triangleright$ $d_{2}=O(k/\epsilon)$

10: Form $\widehat{V}=A_{2}D_{2}(B_{2}D_{2})^{\dagger}\in\mathbb{R}^{n\times
k}$.

11: Form
$B_{3}=\widehat{U}^{\top}\odot\widehat{V}^{\top}\in\mathbb{R}^{k\times
n^{2}}$.

12:
$D_{3}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{3}^{\top},B_{3}^{\top},n^{2},n,k,\epsilon$).
$\triangleright$ Algorithm 7

13: $d_{3}$ denote the number of nonzero entries in $D_{3}$. $\triangleright$
$d_{3}=O(k/\epsilon)$

14: $C\leftarrow A_{1}D_{1}$, $R\leftarrow A_{2}D_{2}$, $T\leftarrow
A_{3}D_{3}$.

15:
$U\leftarrow\sum_{i=1}^{k}((B_{1}D_{1})^{\dagger})_{i}\otimes((B_{2}D_{2})^{\dagger})_{i}\otimes((B_{3}D_{3})^{\dagger})_{i}$.

16: return $C$, $R$, $T$ and $U$.

17:end procedure

######  Theorem C.41.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, let $k\geq
1$, and let $U_{B},V_{B},W_{B}\in\mathbb{R}^{n\times k}$ denote a rank-$k$,
$\alpha$-approximation to $A$. Then there exists an algorithm which takes
$O(\operatorname{nnz}(A)\log n+n^{2}\operatorname{poly}(\log n,k,1/\epsilon))$
time and outputs three matrices: $C\in\mathbb{R}^{n\times c}$ with columns
from $A$, $R\in\mathbb{R}^{n\times r}$ with rows from $A$,
$T\in\mathbb{R}^{n\times t}$ with tubes from $A$, and a tensor
$U\in\mathbb{R}^{c\times r\times t}$ with $\operatorname{rank}(U)=k$ such that
$c=r=t=O(k/\epsilon)$, and

$\displaystyle\left\|\sum_{i=1}^{c}\sum_{j=1}^{r}\sum_{l=1}^{t}U_{i,j,l}\cdot
C_{i}\otimes R_{j}\otimes
T_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\alpha\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}$

holds with probability $9/10$.

###### Proof.

The proof is almost the same as the proof of Theorem C.40. The only difference
is that instead of using Theorem C.38, we use Theorem C.14. ∎

###  C.8 Face-based selection and decomposition

Previously we provided column-based tensor CURT algorithms, which are
algorithms that can select a subset of columns from each of the three
dimensions. Here we provide two face-based tensor CURT decomposition
algorithms. The first algorithm runs in polynomial time and is a bicriteria
algorithm (the number of samples is $\operatorname{poly}(k/\epsilon)$). The
second algorithm needs to start with a rank-$k$ $(1+O(\epsilon))$-approximate
solution, which we then show how to combine with our previous algorithm. Both
of our algorithms are able to select a subset of column-row faces, a subset of
row-tube faces and a subset of column-tube faces. The second algorithm is able
to output $U$, but the first algorithm is not.

####  C.8.1 Column-row, column-tube, row-tube face subset selection

Algorithm 13 Frobenius Norm Tensor Column-row, Row-tube and Tube-column Face
Subset Selection

1:procedure FFaceCRTSelection($A,n,k,\epsilon$) $\triangleright$ Theorem C.42

2: $s_{1}\leftarrow s_{2}\leftarrow O(k/\epsilon)$.

3: Choose a Gaussian matrix $S_{1}$ with $s_{1}$ columns. $\triangleright$
Definition B.18

4: Choose a Gaussian matrix $S_{2}$ with $s_{2}$ columns. $\triangleright$
Definition B.18

5: Form matrix $V_{3}$ by setting the $(i,j)$-th column to be
$(A_{2}S_{2})_{j}$.

6:
$D_{3}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{2}$,$V_{3}$,$n$,$n^{2}$,$s_{1}s_{2}$,$\epsilon$).
$\triangleright$ Algorithm 7

7: Let $d_{3}$ denote the number of nonzero entries in $D_{3}$.
$\triangleright$ $d_{3}=O(s_{1}s_{2}/\epsilon)$

8: Form matrix $U_{2}$ by setting the $(i,j)$-th column to be
$(A_{1}S_{1})_{i}$.

9:
$D_{2}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{1}$,$U_{2}$,$n$,$n^{2}$,$s_{1}s_{2}$,$\epsilon$).

10: Let $d_{2}$ denote the number of nonzero entries in $D_{2}$.
$\triangleright$ $d_{2}=O(s_{1}s_{2}/\epsilon)$

11: Form matrix $W_{1}$ by setting the $(i,j)$-th column to be
$(A(I,D_{3},I)_{3})_{j}$.

12:
$D_{1}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{3}$,$W_{1}$,$n$,$n^{2}$,$s_{1}s_{2}$,$\epsilon$).

13: Let $d_{1}$ denote the number of nonzero entries in $D_{1}$.
$\triangleright$ $d_{1}=O(s_{1}s_{2}/\epsilon)$

14: $T\leftarrow A(I,I,D_{1})$, $C\leftarrow A(D_{2},I,I)$, and $R\leftarrow
A(I,D_{3},I)$.

15: return $C$, $R$ and $T$.

16:end procedure

######  Theorem C.42.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any
$k\geq 1$, there exists an algorithm which takes $O(\operatorname{nnz}(A))\log
n+n^{2}\operatorname{poly}(\log n,k,1/\epsilon)$ time and outputs three
tensors : a subset $C\in\mathbb{R}^{c\times n\times n}$ of row-tube faces of
$A$, a subset $R\in\mathbb{R}^{n\times r\times n}$ of column-tube faces of
$A$, and a subset $T\in\mathbb{R}^{n\times n\times t}$ of column-row faces of
$A$, where $c=r=t=\operatorname{poly}(k,1/\epsilon)$, and for which there
exists a tensor $U\in\mathbb{R}^{tn\times cn\times rn}$ for which

$\displaystyle\|U(T_{1},C_{2},R_{3})-A\|_{F}^{2}\leq(1+\epsilon)\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2},$

or equivalently,

$\displaystyle\left\|\sum_{i=1}^{tn}\sum_{j=1}^{cn}\sum_{l=1}^{rn}U_{i,j,l}\cdot(T_{1})_{i}\otimes(C_{2})_{j}\otimes(R_{3})_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}.$

###### Proof.

We fix $V^{*}\in\mathbb{R}^{n\times k}$ and $W^{*}\in\mathbb{R}^{n\times k}$.
We define $Z_{1}\in\mathbb{R}^{k\times n^{2}}$ where the $i$-th row of $Z_{1}$
is the vector $V_{i}\otimes W_{i}$. Choose a sketching (Gaussian) matrix
$S_{1}\in\mathbb{R}^{n^{2}\times s_{1}}$ (Definition B.18), and let
$\widehat{U}=A_{1}S_{1}(Z_{1}S_{1})^{\dagger}\in\mathbb{R}^{n\times k}$.
Following a similar argument as in the previous theorem, we have

$\displaystyle\|\widehat{U}Z_{1}-A_{1}\|_{F}^{2}\leq(1+\epsilon)\operatorname{OPT}.$

We fix $\widehat{U}$ and $W^{*}$. We define $Z_{2}\in\mathbb{R}^{k\times
n^{2}}$ where the $i$-th row of $Z_{2}$ is the vector $\widehat{U}_{i}\otimes
W^{*}_{i}$. Choose a sketching (Gaussian) matrix
$S_{2}\in\mathbb{R}^{n^{2}\times s_{2}}$ (Definition B.18), and let
$\widehat{V}=A_{2}S_{2}(Z_{2}S_{2})^{\dagger}\in\mathbb{R}^{n\times k}$.
Following a similar argument as in the previous theorem, we have

$\displaystyle\|\widehat{V}Z_{2}-A_{2}\|_{F}^{2}\leq(1+\epsilon)^{2}\operatorname{OPT}.$

We fix $\widehat{U}$ and $\widehat{V}$. Note that
$\widehat{U}=A_{1}S_{1}(Z_{1}S_{1})^{\dagger}$ and
$\widehat{V}=A_{2}S_{2}(Z_{2}S_{2})^{\dagger}$. We define
$Z_{3}\in\mathbb{R}^{k\times n^{2}}$ such that the $i$-th row of $Z_{3}$ is
the vector $\widehat{U}_{i}\otimes\widehat{V}_{i}$. Let $z_{3}=s_{1}\cdot
s_{2}$. We define $Z^{\prime}_{3}\in\mathbb{R}^{z_{3}\times n^{2}}$ such that,
$\forall i\in[s_{1}],\forall j\in[s_{2}]$, the $i+(j-1)s_{1}$-th row of
$Z^{\prime}_{3}$ is the vector $(A_{1}S_{1})_{i}\otimes(A_{2}S_{2})_{j}$.

We define $U_{3}\in\mathbb{R}^{n\times z_{3}}$ to be the matrix where the
$i+(j-1)s_{1}$-th column is $(A_{1}S_{1})_{i}$ and
$V_{3}\in\mathbb{R}^{n\times z_{3}}$ to be the matrix where the
$i+(j-1)s_{1}$-th column is $(A_{2}S_{2})_{j}$. Then
$Z_{3}^{\prime}=(U_{3}^{\top}\odot V_{3}^{\top})$.

We first have,

$\displaystyle\min_{W\in\mathbb{R}^{n\times k},X\in\mathbb{R}^{k\times
z_{3}}}\|WXZ_{3}^{\prime}-A_{3}\|_{F}^{2}\leq\min_{W\in\mathbb{R}^{n\times
k}}\|WZ_{3}-A_{3}\|_{F}^{2}\leq(1+\epsilon)^{2}\operatorname{OPT}.$

Now consider the following objective function,

$\displaystyle\min_{W\in\mathbb{R}^{n\times z_{3}}}\|V_{3}\cdot(W^{\top}\odot
U_{3}^{\top})-A_{2}\|_{F}^{2}.$

Let $D_{3}$ denote a sampling and rescaling diagonal matrix according to
$V_{1}\in\mathbb{R}^{n\times z_{3}}$, let $d_{3}$ denote the number of nonzero
entries of $D_{3}$. Then we have

$\displaystyle\leavevmode\nobreak\ \min_{W\in\mathbb{R}^{n\times
z_{3}}}\|D_{3}V_{3}\cdot(W^{\top}\odot U_{3}^{\top})-D_{3}A_{2}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{W\in\mathbb{R}^{n\times z_{3}}}\|U_{3}\otimes(D_{3}V_{3})\otimes
W-A(I,D_{3},I)\|_{F}^{2}$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{W\in\mathbb{R}^{n\times
z_{3}}}\|W\cdot(U_{3}^{\top}\odot(D_{3}V_{3})^{\top})-(A(I,D_{3},I))_{3}\|_{F}^{2},$

where the first equality follows by retensorizing the objective function, and
the second equality follows by flattening the tensor along the third
dimension.

Let $\overline{Z}_{3}$ denote
$(U_{3}^{\top}\odot(D_{3}V_{3})^{\top})\in\mathbb{R}^{z_{3}\times nd_{3}}$ and
$W^{\prime}=(A(I,D_{3},I))_{3}\in\mathbb{R}^{n\times nd_{3}}$. Using Theorem
C.14, we can find a diagonal matrix $D_{3}\in\mathbb{R}^{n^{2}\times n^{2}}$
with $d_{3}=O(z_{3}/\epsilon)=O(k^{2}/\epsilon^{3})$ nonzero entries such that

$\displaystyle\|U_{3}\otimes
V_{3}\otimes(W^{\prime}Z_{3}^{\dagger})-A\|_{F}^{2}\leq(1+\epsilon)^{3}\operatorname{OPT}.$

We define $U_{2}=U_{3}\in\mathbb{R}^{n\times z_{2}}$ with $z_{2}=z_{3}$. We
define $W_{2}=W^{\prime}\overline{Z}_{3}^{\dagger}\in\mathbb{R}^{n\times
z_{2}}$ with $z_{2}=z_{3}$. We consider,

$\displaystyle\min_{V\in\mathbb{R}^{n\times z_{2}}}\|U_{2}\cdot(V^{\top}\odot
W_{2}^{\top})-A_{1}\|_{F}^{2}.$

Let $D_{2}$ denote a sampling and rescaling matrix according to $U_{2}$, and
let $d_{2}$ denote the number of nonzero entries of $D_{2}$. Then, we have

$\displaystyle\leavevmode\nobreak\ \min_{V\in\mathbb{R}^{n\times
z_{2}}}\|D_{2}U_{2}\cdot(V^{\top}\odot W_{2}^{\top})-D_{2}A_{1}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{V\in\mathbb{R}^{n\times z_{2}}}\|D_{2}U_{2}\otimes V\otimes
W_{2}-A(D_{2},I,I)\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \min_{V\in\mathbb{R}^{n\times
z_{2}}}\|V\cdot(W_{2}^{\top}\odot(D_{2}U_{2})^{\top})-(A(D_{2},I,I))_{2}\|_{F}^{2},$

where the first equality follows by retensorizing the objective function, and
the second equality follows by flattening the tensor along the second
dimension.

Let $\overline{Z}_{2}$ denote
$(W_{2}^{\top}\odot(D_{2}U_{2})^{\top})\in\mathbb{R}^{z_{2}\times nd_{2}}$ and
$V^{\prime}=(A(D_{2},I,I))_{2}\in\mathbb{R}^{n\times nd_{2}}$. Using Theorem
C.14, we can find a diagonal matrix $D_{2}\in\mathbb{R}^{n^{2}\times n^{2}}$
with $d_{2}=O(z_{2}/\epsilon)$ nonzero entries such that

$\displaystyle\|U_{2}\otimes(V^{\prime}\overline{Z}_{2}^{\dagger})\otimes
W_{2}-A\|_{F}^{2}\leq(1+\epsilon)^{4}\operatorname{OPT}.$

We define $W_{1}=W_{2}\in\mathbb{R}^{n\times z_{1}}$ with $z_{1}=z_{2}$, and
define $V_{1}=(V^{\prime}\overline{Z}_{2}^{\dagger})\in\mathbb{R}^{n\times
z_{1}}$ with $z_{1}=z_{2}$.

Let $D_{1}$ denote a sampling and rescaling matrix according to $W_{1}$, and
let $d_{1}$ denote the number of nonzero entries of $D_{1}$. Then we have

$\displaystyle\leavevmode\nobreak\ \min_{U\in\mathbb{R}^{n\times
z_{1}}}\|D_{1}W_{1}\cdot(U^{\top}\odot V_{1}^{\top})-D_{1}A_{3}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{U\in\mathbb{R}^{n\times z_{1}}}\|U\otimes
V_{1}\otimes(D_{1}W_{1})-A(I,I,D_{1})\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \min_{U\in\mathbb{R}^{n\times
z_{1}}}\|U\cdot(V_{1}^{\top}\odot(D_{1}W_{1})^{\top})-A(I,I,D_{1})_{1}\|_{F}^{2}$

where the first equality follows by unflattening the objective function, and
second equality follows by flattening the tensor along the first dimension.

Let $\overline{Z}_{1}$ denote
$(V_{1}^{\top}\odot(D_{1}W_{1})^{\top})\in\mathbb{R}^{z_{1}\times nd_{1}}$,
and $U^{\prime}=A(I,I,D_{1})_{1}\in\mathbb{R}^{n\times nd_{1}}$. Using Theorem
C.14, we can find a diagonal matrix $D_{1}\in\mathbb{R}^{n^{2}\times n^{2}}$
with $d_{1}=O(z_{1}/\epsilon)$ nonzero entries such that

$\displaystyle\|(U^{\prime}\overline{Z}_{1}^{\dagger})\otimes(V_{1})\otimes
W_{1}-A\|_{F}^{2}\leq(1+\epsilon)^{5}\operatorname{OPT},$

which means,

$\displaystyle\|(U^{\prime}\overline{Z}_{1}^{\dagger})\otimes(V^{\prime}\overline{Z}_{2}^{\dagger})\otimes(W^{\prime}\overline{Z}_{3}^{\dagger})-A\|_{F}^{2}\leq(1+\epsilon)^{5}\operatorname{OPT}.$

Putting $U^{\prime},V^{\prime},W^{\prime}$ together completes the proof. ∎

######  Corollary C.43.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, for any
$k\geq 1$, there exists an algorithm which takes
$O(\operatorname{nnz}(A))+n^{2}\operatorname{poly}(k,1/\epsilon)$ time and
outputs three tensors : a subset $C\in\mathbb{R}^{c\times n\times n}$ of row-
tube faces of $A$, a subset $R\in\mathbb{R}^{n\times r\times n}$ of column-
tube faces of $A$, and a subset $T\in\mathbb{R}^{n\times n\times t}$ of
column-row faces of $A$, where $c=r=t=\operatorname{poly}(k,1/\epsilon)$, so
that there exists a tensor $U\in\mathbb{R}^{tn\times cn\times rn}$ for which

$\displaystyle\|U(T_{1},C_{2},R_{3})-A\|_{F}^{2}\leq(1+\epsilon)\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2},$

or equivalently,

$\displaystyle\left\|\sum_{i=1}^{tn}\sum_{j=1}^{cn}\sum_{l=1}^{rn}U_{i,j,l}\cdot(T_{1})_{i}\otimes(C_{2})_{j}\otimes(R_{3})_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}$

###### Proof.

If we allow a $\operatorname{poly}(k/\epsilon)$ factor increase in running
time and a $\operatorname{poly}(k/\epsilon)$ factor increase in the number of
faces selected, then instead of using generalized row subset selection, which
has running time depending on $\log n$, we can use the technique in Section
C.11 to avoid the $\log n$ factor. ∎

####  C.8.2 CURT decomposition

Algorithm 14 Frobenius Norm (Face-based) CURT Decomposition Algorithm, Optimal
Sample Complexity

1:procedure FFaceCURTDecomposition($A,U_{B},V_{B},W_{B},n,k$) $\triangleright$
Theorem C.44

2:
$D_{1}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{3},W_{B},n,n^{2},k,\epsilon$).
$\triangleright$ Algorithm 7, the number of nonzero entries is
$d_{1}=O(k/\epsilon)$

3: Form $Z_{1}=V_{B}^{\top}\odot(D_{1}W_{B})^{\top}$.

4: Form $\widehat{U}=(A(I,I,D_{1}))_{1}Z_{1}^{\dagger}\in\mathbb{R}^{n\times
k}$.

5:
$D_{2}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{1},\widehat{U},n,n^{2},k,\epsilon$).
$\triangleright$ The number of nonzero entries is $d_{2}=O(k/\epsilon)$

6: Form $Z_{2}=(W_{B}^{\top}\odot(D_{2}\widehat{U}))$.

7: Form $\widehat{V}=(A(D_{2},I,I))_{2}Z_{2}^{\dagger}\in\mathbb{R}^{n\times
k}$.

8:
$D_{3}\leftarrow$GeneralizedMatrixRowSubsetSelection($A_{2},\widehat{V},n,n^{2},k,\epsilon$).
$\triangleright$ The number of nonzero entries is $d_{3}=O(k/\epsilon)$

9: Form $Z_{3}=\widehat{U}^{\top}\odot(D_{3}\widehat{V})^{\top}$.

10: Form
$\widehat{W}=(A(I,D_{3},I))_{3}(Z_{3})^{\dagger}\in\mathbb{R}^{n\times k}$.

11: $T\leftarrow A(I,I,D_{1})$, $C\leftarrow A(D_{2},I,I)$, $R\leftarrow
A(I,D_{3},I)$.

12:
$U\leftarrow\sum_{i=1}^{k}((Z_{1})^{\dagger})_{i}\otimes((Z_{2})^{\dagger})_{i}\otimes((Z_{3})^{\dagger})_{i}$.

13: return $C$, $R$, $T$ and $U$.

14:end procedure

######  Theorem C.44.

Given a $3$rd order tensor $A\in\mathbb{R}^{n\times n\times n}$, let $k\geq
1$, and let $U_{B},V_{B},W_{B}\in\mathbb{R}^{n\times k}$ denote a rank-$k$,
$\alpha$-approximation to $A$. Then there exists an algorithm which takes
$O(\operatorname{nnz}(A))\log n+n^{2}\operatorname{poly}(\log n,k,1/\epsilon)$
time and outputs three tensors: $C\in\mathbb{R}^{c\times n\times n}$ with row-
tube faces from $A$, $R\in\mathbb{R}^{n\times r\times n}$ with colum-tube
faces from $A$, $T\in\mathbb{R}^{n\times n\times t}$ with column-row faces
from $A$, and a (factorization of a) tensor $U\in\mathbb{R}^{tn\times cn\times
rn}$ with $\operatorname{rank}(U)=k$ for which $c=r=t=O(k/\epsilon)$ and

$\displaystyle\|U(T_{1},C_{2},R_{3})-A\|_{F}^{2}\leq(1+\epsilon)\alpha\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2},$

or equivalently,

$\displaystyle\left\|\sum_{i=1}^{tn}\sum_{j=1}^{cn}\sum_{l=1}^{rn}U_{i,j,l}\cdot(T_{1})_{i}\otimes(C_{2})_{j}\otimes(R_{3})_{l}-A\right\|_{F}^{2}\leq(1+\epsilon)\alpha\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}$

holds with probability $9/10$.

###### Proof.

We already have three matrices $U_{B}\in\mathbb{R}^{n\times k}$,
$V_{B}\in\mathbb{R}^{n\times k}$ and $W_{B}\in\mathbb{R}^{n\times k}$ and
these three matrices provide a $\operatorname{rank}$-$k$,
$\alpha$-approximation to $A$, i.e.,

$\displaystyle\|U_{B}\otimes V_{B}\otimes
W_{B}-A\|_{F}^{2}\leq\alpha\underbrace{\min_{\operatorname{rank}-k\leavevmode\nobreak\
A^{\prime}}\|A^{\prime}-A\|_{F}^{2}}_{\operatorname{OPT}}.$

We can consider the following problem,

$\displaystyle\min_{U\in\mathbb{R}^{n\times k}}\|W_{B}\cdot(U^{\top}\odot
V_{B}^{\top})-A_{3}\|_{F}^{2}.$

Let $D_{1}$ denote a sampling and rescaling diagonal matrix according to
$W_{B}$, and let $d_{1}$ denote the number of nonzero entries of $D_{1}$. Then
we have

$\displaystyle\min_{U\in\mathbb{R}^{n\times
k}}\|(D_{1}W_{B})\cdot(U^{\top}\odot V_{B}^{\top})-D_{1}A_{3}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\min_{U\in\mathbb{R}^{n\times k}}\|U\otimes
V_{B}\otimes D_{1}W_{B}-A(I,I,D_{1})\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\min_{U\in\mathbb{R}^{n\times
k}}\|U\cdot(V_{B}^{\top}\odot(D_{1}W_{B})^{\top})-(A(I,I,D_{1}))_{1}\|_{F}^{2},$

where the first equality follows by retensorizing the objective function, and
the second equality follows by flattening the tensor along the first
dimension. Let $Z_{1}$ denote
$V_{B}^{\top}\odot(D_{1}W_{B})^{\top}\in\mathbb{R}^{k\times nd_{1}}$, and
define $\widehat{U}=(A(I,I,D_{1}))_{1}Z_{1}^{\dagger}\in\mathbb{R}^{n\times
k}$. Then we have

$\displaystyle\|\widehat{U}\otimes V_{B}\otimes
W_{B}-A\|_{F}^{2}\leq(1+\epsilon)\alpha\operatorname{OPT}.$

In the second step, we fix $\widehat{U}$ and $W_{B}$, and consider the
following objective function,

$\displaystyle\min_{V\in\mathbb{R}^{n\times
k}}\|\widehat{U}\cdot(V^{\top}\odot W_{B})-A_{1}\|_{F}^{2}.$

Let $D_{2}$ denote a sampling and rescaling matrix according to $\widehat{U}$,
and let $d_{2}$ denote the number of nonzero entries of $D_{2}$. Then we have,

$\displaystyle\leavevmode\nobreak\ \min_{V\in\mathbb{R}^{n\times
k}}\|(D_{2}\widehat{U})\cdot(V^{\top}\odot W_{B}^{\top})-D_{2}A_{1}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{V\in\mathbb{R}^{n\times k}}\|(D_{2}\widehat{U})\otimes V\otimes
W_{B}-A(D_{2},I,I)\|_{F}^{2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \min_{V\in\mathbb{R}^{n\times
k}}\|V\cdot(W_{B}^{\top}\odot(D_{2}\widehat{U})^{\top})-(A(D_{2},I,I))_{2}\|_{F}^{2},$

where the first equality follows by unflattening the objective function, and
the second equality follows by flattening the tensor along the second
dimension. Let $Z_{2}$ denote
$(W_{B}^{\top}\odot(D_{2}\widehat{U})^{\top})\in\mathbb{R}^{k\times nd_{2}}$,
and define
$\widehat{V}=(A(D_{2},I,I))_{2}(Z_{2})^{\dagger}\in\mathbb{R}^{n\times k}$.
Then we have,

$\displaystyle\|\widehat{U}\otimes\widehat{V}\otimes
W_{B}-A\|_{F}^{2}\leq(1+\epsilon)^{2}\alpha\operatorname{OPT}.$

In the third step, we fix $\widehat{U}$ and $\widehat{V}$, and consider the
following objective function,

$\displaystyle\min_{W\in\mathbb{R}^{n\times
k}}\|\widehat{V}\cdot(W\odot\widehat{U})-A_{2}\|_{F}^{2}.$

Let $D_{3}$ denote a sampling and rescaling matrix according to $\widehat{V}$,
and let $d_{3}$ denote the number of nonzero entries of $D_{3}$. Then we have,

$\displaystyle\leavevmode\nobreak\ \min_{W\in\mathbb{R}^{n\times
k}}\|(D_{3}\widehat{V})\cdot(W^{\top}\odot\widehat{U}^{\top})-D_{3}A_{2}\|_{F}^{2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{W\in\mathbb{R}^{n\times k}}\|\widehat{U}\otimes(D_{3}\widehat{V})\otimes
W-A(I,D_{3},I)\|_{F}^{2}$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\min_{W\in\mathbb{R}^{n\times
k}}\|W\cdot(\widehat{U}^{\top}\odot(D_{3}\widehat{V})^{\top})-(A(I,D_{3},I))_{3}\|_{F}^{2},$

where the first equality follows by retensorizing the objective function, and
the second equality follows by flattening the tensor along the third
dimension. Let $Z_{3}$ denote
$(\widehat{U}^{\top}\odot(D_{3}\widehat{V})^{\top})\in\mathbb{R}^{k\times
nd_{3}}$, and define $\widehat{W}=(A(I,D_{3},I))_{3}(Z_{3})^{\dagger}$.
Putting it all together, we have,

$\displaystyle\|\widehat{U}\otimes\widehat{V}\otimes\widehat{W}-A\|_{F}^{2}\leq(1+\epsilon)^{3}\alpha\operatorname{OPT}.$

This implies

$\displaystyle\|(A(I,I,D_{1}))_{1}Z_{1}^{\dagger}\otimes(A(D_{2},I,I))_{2}Z_{2}^{\dagger}\otimes(A(I,D_{3},I))_{3}Z_{3}^{\dagger}-A\|_{F}^{2}\leq(1+\epsilon)^{3}\alpha\operatorname{OPT}.$

∎

###  C.9 Solving small problems

######  Theorem C.45.

Let $\max_{i}\\{t_{i},d_{i}\\}\leq n$. Given a $t_{1}\times t_{2}\times t_{3}$
tensor $A$ and three matrices: a $t_{1}\times d_{1}$ matrix $T_{1}$, a
$t_{2}\times d_{2}$ matrix $T_{2}$, and a $t_{3}\times d_{3}$ matrix $T_{3}$,
if for any $\delta>0$ there exists a solution to

$\displaystyle\min_{X_{1},X_{2},X_{3}}\left\|\sum_{i=1}^{k}(T_{1}X_{1})_{i}\otimes(T_{2}X_{2})_{i}\otimes(T_{3}X_{3})_{i}-A\right\|_{F}^{2}:=\operatorname{OPT},$

and each entry of $X_{i}$ can be expressed using $O(n^{\delta})$ bits, then
there exists an algorithm that takes $n^{O(\delta)}\cdot
2^{O(d_{1}k+d_{2}k+d_{3}k)}$ time and outputs three matrices:
$\widehat{X}_{1}$, $\widehat{X}_{2}$, and $\widehat{X}_{3}$ such that
$\|(T_{1}\widehat{X}_{1})\otimes(T_{2}\widehat{X}_{2})\otimes(T_{3}\widehat{X}_{3})-A\|_{F}^{2}=\operatorname{OPT}$.

###### Proof.

For each $i\in[3]$, we can create $t_{i}\times d_{i}$ variables to represent
matrix $X_{i}$. Let $x$ denote this list of variables. Let $B$ denote tensor
$\sum_{i=1}^{k}(T_{1}X_{1})_{i}\otimes(T_{2}X_{2})_{i}\otimes(T_{3}X_{3})_{i}$
and let $B_{i,j,l}(x)$ denote an entry of tensor $B$ (which can be thought of
as a polynomial written in terms of $x$). Then we can write the following
objective function,

$\displaystyle\min_{x}\sum_{i=1}^{t_{1}}\sum_{j=1}^{t_{2}}\sum_{l=1}^{t_{3}}(B_{i,j,l}(x)-A_{i,j,l})^{2}.$

We slightly modify the above objective function to obtain a new objective
function,

$\displaystyle\min_{x,\sigma}$ $\displaystyle\leavevmode\nobreak\
\sum_{i=1}^{t_{1}}\sum_{j=1}^{t_{2}}\sum_{l=1}^{t_{3}}(B_{i,j,l}(x)-A_{i,j,l})^{2},$
s.t. $\displaystyle\leavevmode\nobreak\ \|x\|_{2}^{2}\leq 2^{O(n^{\delta})},$

where the last constraint is unharmful, because there exists a solution that
can be written using $O(n^{\delta})$ bits. Note that the number of inequality
constraints in the above system is $O(1)$, the degree is $O(1)$, and the
number of variables is $v=(d_{1}k+d_{2}k+d_{3}k)$. Thus by Theorem B.11, the
minimum nonzero cost is at least

$\displaystyle(2^{O(n^{\delta})})^{-2^{{O}(v)}}.$

It is clear that the upper bound on the cost is at most $2^{O(n^{\delta})}$.
Thus the number of binary search steps is at most
$\log(2^{O(n^{\delta})})2^{{O}(v)}$. In each step of the binary search, we
need to choose a cost $C$ between the lower bound and the upper bound, and
write down the polynomial system,

$\displaystyle\leavevmode\nobreak\
\sum_{i=1}^{t_{1}}\sum_{j=1}^{t_{2}}\sum_{l=1}^{t_{3}}(B_{i,j,l}(x)-A_{i,j,l})^{2}\leq
C,$ $\displaystyle\leavevmode\nobreak\ \|x\|_{2}^{2}\leq 2^{O(n^{\delta})}.$

Using Theorem B.10, we can determine if there exists a solution to the above
polynomial system. Since the number of variables is $v$, and the degree is
$O(1)$, the number of inequality constraints is $O(1)$. Thus, the running time
is

$\displaystyle\operatorname{poly}(\text{bitsize})\cdot(\\#\operatorname{constraints}\cdot\operatorname{degree})^{\\#\operatorname{variables}}=n^{O(\delta)}2^{O(v)}.$

∎

###  C.10 Extension to general $q$-th order tensors

This section provides the details for our extensions from $3$rd order tensors
to general $q$-th order tensors. In most practical applications, the order $q$
is a constant. Thus, to simplify the analysis, we use $O_{q}(\cdot)$ to hide
dependencies on $q$.

####  C.10.1 Fast sampling of columns according to leverage scores, implicitly

This section explains an algorithm that is able to sample from the leverage
scores from the $\odot$ product of $q$ matrices $U_{1},U_{2},\cdots,U_{q}$
without explicitly writing down $U_{1}\odot U_{2}\odot\cdots U_{q}$. To build
this algorithm we combine TensorSketch, some ideas from [DMIMW12], and some
techniques from [AKO11, MW10]. Finally, we improve the running time for
sampling columns according to the leverage scores from
$\operatorname{poly}(n)$ to $\widetilde{O}(n)$. Given $q$ matrices
$U_{1},U_{2},\cdots,U_{q}$, with each such matrix $U_{i}$ having size $k\times
n_{i}$, we define $A\in\mathbb{R}^{k\times\prod_{i=1}^{q}n_{i}}$ to be the
matrix where the $i$-th row of $A$ is the vectorization of $U_{1}^{i}\otimes
U_{2}^{i}\otimes\cdots\otimes U_{q}^{i}$, $\forall i\in[k]$. Naïvely, in order
to sample $\operatorname{poly}(k,1/\epsilon)$ rows from $A$ according to the
leverage scores, we need to write down $\prod_{i=1}^{q}n_{i}$ leverage scores.
This approach will take at least $\prod_{i=1}^{q}n_{i}$ running time. In the
remainder of this section, we will explain how to do it in