# Learning Sparsity and Randomness for Data-driven Low Rank Approximation

Tiejin Chen∗
Department of Statistics
University of Michigan
<EMAIL_ADDRESS>
&Yicheng Tao∗
Department of Computer Science
University of Michigan
<EMAIL_ADDRESS>

###### Abstract

Learning-based low rank approximation algorithms can significantly improve the
performance of randomized low rank approximation with sketch matrix. With the
learned value and fixed non-zero positions for sketch matrices from learning-
based algorithms, these matrices can reduce the test error of low rank
approximation significantly. However, there is still no good method to learn
non-zero positions as well as overcome the out-of-distribution performance
loss.

In this work, we introduce two new methods Learning Sparsity and Learning
Randomness which try to learn a better sparsity patterns and add randomness to
the value of sketch matrix. These two methods can be applied with any
learning-based algorithms which use sketch matrix directly. Our experiments
show that these two methods can improve the performance of previous learning-
based algorithm for both test error and out-of-distribution test error without
adding too much complexity.

**footnotetext: Equal contribution.

##  1 Introduction

The problem of low-rank approximation (LRA) of matrices is that given a matrix
$A\in\mathbb{R}^{n\times d}$ with $n\geq d$ and an integer $0<k\ll d$, find a
rank-$k$ matrix $A^{\prime}$ that minimizes the approximation error
$\left\lVert A-A^{\prime}\right\rVert_{F}^{2}$ . Truncated Singular value
decomposition (SVD) can solve this problem optimally in polynomial time, which
is nevertheless too slow for large matrices. To handle this, fast approximate
LRA algorithms based on sketching have been proposed [14, 6, 3, 15], which
basically construct a sketch matrix $S\in\mathbb{R}^{m\times n}$ and compute
LRA from a much smaller matrix $SA$ to get a save of both time and storage
space. However, these algorithms suffer from a larger approximation error than
the optimal LRA, which may affect downstream tasks based on LDA. So the
challenge here is to design a sketch matrix that can also produce a good
approximation quality.

Data-driven algorithms have been recently developed to handle this challenge
[12, 1, 17, 13]. They leverage past data that are related to future input as
training samples to learn sketching matrices, which result in significant
improvement over non-data-driven baselines. However, as mentioned in [13], the
sparsity pattern of the sketching matrix has not been considered to be learned
during training, which may additionally improve performance. All the
algorithms in the this area consider to learn fixed value for sketch matrix
which may cause the problem of over-fitting, and result in poor performance in
out-of-distribution situation. Indyk et al.[12] provide a proof of worst case
bound which concatenates learned sketch matrix with data-oblivious random
matrix vertically. Nonetheless, such kind of method will increase sketch size
and get a reduce of efficiency.

In this work, we propose two different methods to overcome the disadvantages
we mentioned above. The first method Learning Sparsity will learn the sparsity
patterns given target total number of non-zero positions as well as value of
sketch matrix $S$. We will use an extra position matrix $D$ to learn the
position, and use $S^{D}=S\odot D$ as final sketch matrix where $\odot$
represents Hardmard product. The second method Learning Randomness tends to
learn a Gaussian distribution for every non-zero position of sketch matrix
instead of a fixed learned value. We will learn the mean and variance of the
distribution. It is easy to notice that our two methods can be easily applied
to any learning-based LRA algorithms require sketch matrix $S$ directly.
Finally, we combine our two methods to get a new algorithm based on IVY [12].

Overall our main contribution can be summarized as follows:

  * •

We propose two flexible methods Learning Sparsity and Learning Randomness
which can be easily applied to any learning-based LRA method requires sketch
matrix $S$ directly.

  * •

We combine these two methods to get a new algorithm based on IVY. Our new
algorithm can improve the performance of IVY significantly.

  * •

We conduct experiments to show our method can actually improve the performance
of learning-based algorithms.

##  2 Related Work

Low rank approximation with sketch matrix sample from different random
distributions [9, 5, 20, 11] has been applied wildly but suffers from high
approximation error compared to optimal solution. Learning-based LRA is first
introduced in [12] and it proposes IVY, a learning-based LRA algorithm
directly follows a common streaming algorithm with random sketch matrix named
SCW [4, 19]. IVY can increase the performance of sketch methods significantly.
After that, some work try to increase the performance of IVY. Liu et al. [17]
tends to learn a better sparsity patterns for CountSketch by trying every
position after every iteration during training. FewshotSGD [13] is proposed to
reduce the training time with surrogate loss instead of direct empirical loss
used in IVY. Most recently, Liu et al. [16] come up with a tensor-based
algorithm with tensor decomposition while Sakaue and Oki [18] demonstrate
iterative hard thresholding method to learn better non-zero positions for
sketch matrix.

There are only few work related to the theoretical part of learning-based LRA
especially considering error bound. Bartlett et al. [2] give a generalization
bound for learning-based methods and this result can be improved [18] by PAC-
learning approach [7, 8].

##  3 Preliminaries

Given a data matrix $A\in\mathbb{R}^{n\times d}$ and a sketching matrix
$S\in\mathbb{R}^{m\times n}$, SCW [4] is a widely-used algorithm for sketch-
based LDA [12, 13, 17], which is shown in Algorithm 1. It first computes the
compact SVD of the sketch $SA$ as $U\Sigma V^{\top}$, then computes the best
rank-$k$ approximation of $AV$ denoted as $[AV]_{k}$, and finally outputs
$[AV]_{k}V^{\top}$ as a rank-$k$ approximation of $A$. We remark that SCW is
differential if it uses a differentiable SVD algorithm, which is available in
PyTorch and used by us.

Input: data matrix $A\in\mathbb{R}^{n\times d}$, sketching matrix
$S\in\mathbb{R}^{m\times n}$

1

2$U,\Sigma,V^{\top}\leftarrow\text{CompactSVD}(SA)\quad\triangleright\quad\left\\{r=\operatorname{rank}(SA),U\in\mathbb{R}^{m\times
r},V\in\mathbb{R}^{d\times r}\right\\}$

Output: $[AV]_{k}V^{\top}$

Algorithm 1 SCW

##  4 Method

We propose two new computational components, namely learning-based sparsity
and randomness, for IVY. The learning-based sparsity component provides a
solution to the challenge of learning the sparsity pattern of the sketching
matrix, while the learning-based randomness component embeds learned Gaussian
randomness in the formation of the sketching matrix to facilitate the
generalization ability of IVY. We combine each component as well as both of
them with IVY and get three new algorithms IVY+LS, IVY+LR, and IVY+LS&LR.

###  4.1 Learning-based Sparsity

In the previous works of data-driven LDA, the non-zero values of the sketching
matrix are learned with their positions fixed, which are randomly chosen
before training. This obviously limits the model from learning the best
positions for non-zero entries. We fix this issue by learning a sparse matrix
$D\in\mathbb{R}_{+}^{m\times n}$ that indicates the sparsity pattern of the
sketching matrix $S$. Specifically, we element-wisely multiply $D$ with $S$ to
make $S$ sparse and input the new $S$ to the SCW algorithm. During training,
we use an approximation error loss with an $L_{1}$ regularization on $D$ to
control the sparsity of the learned $D$, and if $D$’s entries are less than a
threshold $\epsilon$, they will be immediately set to zeros. After $S$ meets
our requirement for sparsity, we will stop train $D$ but continue to train
$S$. The algorithm of this component combined with IVY, IVY+LS, is shown in
Algorithm 2.

Input: training set $\\{A_{1},\cdots,A_{N}\\}\subset\mathbb{R}^{n\times d}$,
sketch size $m$, target sparsity $s$, regularization weight $\lambda$,
threshold $\epsilon$, learning rate $\eta$

1

2$S\leftarrow 1^{m\times n}\in\mathbb{R}^{m\times
n}\quad\triangleright\quad\\{\text{sketching matrix}\\}$

3$D\leftarrow 1^{m\times n}\in\mathbb{R}_{+}^{m\times
n}\quad\triangleright\quad\\{\text{sparsity pattern matrix}\\}$

4for  _$i\leftarrow 1,\cdots,N$_ do

5

6 $S\leftarrow D\otimes S$

7 $L\leftarrow\lVert
A_{i}-\text{SCW}(A_{i},S)\rVert_{F}+\lambda\sum_{i=1}^{m}\sum_{j=1}^{n}D_{ij}$

8 $S\leftarrow S-\eta\frac{\partial L}{\partial S}$

9 if  _$\sum_{i=1}^{m}\sum_{j=1}^{n}1_{D_{ij}\neq 0} >ns$_ then

10

11 $D\leftarrow D-\eta\frac{\partial L}{\partial D}$

12 $D[D<\epsilon]\leftarrow 0$

13 end if

14

15 end for

Output: $S$

Algorithm 2 IVY+LS

###  4.2 Learning-based Randomness

IVY directly learns the values of non-zero entries in the sketching matrix
$S$, which may lead to the overfitting problem and make $S$ less
generalizable. To handle this, we introduce learning-based randomness in the
formation of $S$. Specifically, we see every non-zero entry in $S$ sampled
from an independent Gaussian distribution with learnable mean $\mu$ and
variance $\sigma$ and do sampling in every training iteration. By using the
reparameterization trick
$\mathcal{N}(\mu,\sigma^{2})=\sigma\mathcal{N}(0,1)+\mu$, we can easily update
$\mu$ and $\sigma$ during training. Once learned, we can sample every entry of
$S$ using $\mu$ and $\sigma$ for inference. The algorithm of this component
combined with IVY, IVY+LR, is shown in Algorithm 3. The whole algorithm with
both components combined with IVY, IVY+LS&LR, is shown in Algorithm 4.

Input: training set $\\{A_{1},\cdots,A_{N}\\}\subset\mathbb{R}^{n\times d}$,
sketch size $m$, learning rate $\eta$

1

2$S\leftarrow 0^{m\times n}\in\mathbb{R}^{m\times
n}\quad\triangleright\quad\\{\text{sketching matrix}\\}$

3$\mu,\Sigma\leftarrow 0^{m\times n}\in\mathbb{R}^{m\times n},1^{m\times
n}\in\mathbb{R}_{+}^{m\times n}$

4for  _$i\leftarrow 1,\cdots,N$_ do

5

6 sample $Z\sim\mathcal{N}(0,1/4)^{m\times n}$

7 $S\leftarrow Z\otimes\sqrt{\Sigma}+\mu$

8 $L\leftarrow\lVert A_{i}-\text{SCW}(A_{i},S)\rVert_{F}$

9 $\mu,\Sigma\leftarrow\mu-\eta\frac{\partial
L}{\partial\mu},\text{relu}(\Sigma-\eta\frac{\partial L}{\partial\Sigma})$

10 end for

Output: $S$

Algorithm 3 IVY+LR

Input: training set $\\{A_{1},\cdots,A_{N}\\}\subset\mathbb{R}^{n\times d}$,
sketch size $m$, target sparsity $s$, regularization weight $\lambda$,
threshold $\epsilon$, learning rate $\eta$, early stop flag $flag$

1

2$S\leftarrow 0^{m\times n}\in\mathbb{R}^{m\times
n}\quad\triangleright\quad\\{\text{sketching matrix}\\}$

3$D,\mu,\Sigma\leftarrow 1^{m\times n}\in\mathbb{R}_{+}^{m\times n},0^{m\times
n}\in\mathbb{R}^{m\times n},1^{m\times n}\in\mathbb{R}_{+}^{m\times
n}\quad\triangleright\quad\\{D:\text{sparsity pattern matrix}\\}$

4for  _$i\leftarrow 1,\cdots,N$_ do

5

6 sample $Z\sim\mathcal{N}(0,1/4)^{m\times n}$

7 $S\leftarrow D\otimes(Z\otimes\sqrt{\Sigma}+\mu)$

8 $is\\_sparse\leftarrow\sum_{i=1}^{m}\sum_{j=1}^{n}1_{D_{ij}\neq 0}>ns$

9 if  _is_sparse_ then

10 $\lambda\leftarrow 0$

11 end if

12

13 $L\leftarrow\lVert
A_{i}-\text{SCW}(A_{i},S)\rVert_{F}+\lambda\sum_{i=1}^{m}\sum_{j=1}^{n}\lvert
D_{ij}\rvert$

14 $\mu,\Sigma\leftarrow\mu-\eta\frac{\partial
L}{\partial\mu},\text{relu}(\Sigma-\eta\frac{\partial L}{\partial\Sigma})$

15 if  _$is\\_sparse$ & $flag$_ then

16

17 $D\leftarrow D-\eta\frac{\partial L}{\partial D}$

18 $D[D<\epsilon]\leftarrow 0$

19 end if

20

21 end for

Output: $S$

Algorithm 4 IVY+LS&LR

##  5 Experimental results

###  5.1 Setting

It can be noticed that our two methods Learning Sparsity(LS) and Learning
Randomness(LR) can be applied to any learning-based algorithms requires sketch
matrix $S$ directly besides IVY. In other word, our method can be encode to
many algorithms except Butterfly[1] which does not use $S$ directly. However,
as mentioned in the previous section, we will only encode our methods to IVY
to prove the validity.

Dataset We will only use Hyper dataset from [10] as our training set. Hyper
contains hyperspectral images from natural scenes and each matrix has size of
$1024\times 768$. We will randomly sample 500 images as training set while use
100 images as test set. We will also use Fish dataset† from Kaggle as an
unseen test set during whole training process to verify generalization ability
and out-of-distribution performance of our methods. Fish contains the images
from different classes of fish. We randomly choose 100 images with size of
$1024\times 768$ and transform them into grayscale images.

Sparsity and Density Overall, one measure to describe the sparsity of one
matrix is total number of non-zero entries. However, here we will define
density to describe the sparsity of sketch matrix. More specifically, density
means total number of non-zero positions in one column for sketch matrix.
CountSketch is a matrix with density1 and a matrix with density equal to its
sketch size is a dense matrix.

Metric For the test error, we follow the instruction from IVY [12] and
fewshotSGD[13] and use a test error as:

$err_{A}(S)=\lVert A-SCW(A,S)\rVert_{F}^{2}-\lVert A-A_{k}\rVert_{F}^{2}$ (1)

Where $A_{k}$ is the optimal rank k matrix respect to $A$ and it can be
obtained by truncated SVD.

In all experiments, we will use target rank is 10 which is from the original
paper from IVY. As for training strategies, we follow the original strategy
from IVY with learning rate equals to 1 and momentum equals to 1. The
threshold $\epsilon$ is set to 0.5 and parameter $\lambda$ for $L_{1}$
regularization is set to $0.0003$ when density is 1 and $0.0001$ when density
is larger than 1. All the algorithms are implemented by PyTorch 1.13.0.

$\dagger$$\dagger$footnotetext: can be found at
https://www.kaggle.com/datasets/crowww/a-large-scale-fish-
dataset?resource=download

###  5.2 Average Test Error

Firstly, we test all our method as well as IVY on the Hyper with different
$m$. For Learning Sparsity, we can not direct control density of sketch matrix
while we can control the sparsity e.g the total number of non-zero entries in
sketch matrix. We will report the test error for sketch matrix learned by
Learning Sparsity with less sparsity than other methods. We run the each
setting for 5 times and report the average test error. The results are shown
in Table 1.

As we can see, Learning Sparsity can increase the performance of IVY
significantly regardless of $m$. It seems that our method can learn a better
sparse pattern than random choice. On the other hand, Learning Randomness
cannot perform very well when sketch size is small but do improve the
performance when sketch size is more reasonable. After combining Learning
Sparsity and Learning Randomness, the performance will further increase during
most situation. However, with very small sketch size, the disadvantage of
randomness is still dominating.

Table 1: Test error with density1 $k$, $m$ | IVY | IVY+LS | IVY+LR | IVY+LS&LR
---|---|---|---|---
10, 10 | 4.12 | 3.59 | 6.20 | 5.07
10, 20 | 1.14 | 0.95 | 1.55 | 0.86
10, 40 | 0.43 | 0.23 | 0.41 | 0.21
10, 80 | 0.12 | 0.07 | 0.11 | 0.06

Now let us take a look at the influence of different density. In this
experiment, we will fix $k$ and $m$ to be 10 and 40 respectively. Again, for
Learning Sparsity and combined method, once the sketch matrix meets the
requirement of sparsity we will stop train $D$. We will run in total 500
iterations for the experiments. Our result are shown in Table 2. Overall,
Learning Sparsity can actually learn a better sparse patterns while Learning
Randomness can also work well in some situations. We visualize the results in
Figure 1(a) and Figure 1(b).

(a) Density 1

(b) Density 20

Figure 1: Test error according to training time Table 2: Test error with different density based on sketch size 40 Density | IVY | IVY+LS | IVY+LR | IVY+LS&LR
---|---|---|---|---
1 | 0.43 | 0.23 | 0.41 | 0.21
2 | 0.29 | 0.21 | 0.37 | 0.20
5 | 0.23 | 0.17 | 0.23 | 0.19
10 | 0.20 | 0.13 | 0.21 | 0.17
20 | 0.15 | 0.11 | 0.13 | 0.16
40 | 0.10 | - | 0.08 | -

###  5.3 Running Time

The inference time of our methods will be same as IVY or random matrix since
the output of our methods will be no different compared to them. However,
training time is also considered as one metric to measure the algorithms in
the area. Here we report the training time of 500 iterations for all the
methods including IVY with $m=40$, density$=1$ in Table 3. We can find that
our two new methods will not affect too much about training time. We also
visualize our result for $m=40$, density$=1$ and $m=40$, density $=20$ in
Figure 2(a) and Figure 2(b). Noticing that before sketch matrix from Learning
Sparsity meet our sparsity requirement, the comparison between Learning
Sparsity and other method is unfair. That is the reason why the data for
Learning Sparsity (and combined method) is not start from zero. From the
figures we can easily see that our methods are better even considering
training time instead of iteration numbers.

Table 3: Training time of 500 iterations IVY | IVY+LS | IVY+LR | IVY+LS&LR
---|---|---|---
26.78s | 26.99s | 26.84s | 27.01s

(a) Density 1

(b) Density 20

Figure 2: Test error according to training time

###  5.4 Out-of-distribution Test Error

In this section, we will report our results for the unseen test set Fish to
verify the generalization ability for our methods. We will apply the best
sketch matrix which means it can reach the minimum test error on Hyper to Fish
and get the test error for Fish on $m=40$. The results are shown in Table 4.
We also keep track of test error for Hyper and the results are shown in
Figure3(a) and Figure 3(b). The results proves that our methods can improve
out-of-distribution performance. However we need to be careful if we use
combined algorithm. And we can see from figures, out-of-distribution test
error will decrease with increasing of training time which prove that
learning-based algorithms have generalization ability.

Table 4: Fish Test error ($\times 10^{3}$) with different density based on sketch size 40 Density | IVY | IVY+LS | IVY+LR | IVY+LS&LR
---|---|---|---|---
1 | 0.52 | 0.24 | 0.47 | 0.28
2 | 0.37 | 0.24 | 0.45 | 0.26
5 | 0.29 | 0.19 | 0.31 | 0.24
10 | 0.24 | 0.16 | 0.26 | 0.21
20 | 0.18 | 0.13 | 0.16 | 0.18
40 | 0.09 | - | 0.08 | -

(a) Density 1

(b) Density 20

Figure 3: OOD Test error according to training time

##  6 Conclusion and Further Work

In this paper we propose two new methods to learn the sparsity and randomness
for sketch matrix separately and provide a new algorithm combining this two
methods. It is not hard to notice that both two method can be encode to any
learning-based LRA algorithm with explicit use of sketch matrix. Our
experiments shows our methods can improve the performance of original IVY.
However, in our paper we do not give any theoretical result for neither of our
new method, which we leave for further work. Also, combining these two methods
does not necessarily lead to better performance than each of them, especially
when the sketch size is large, which is another interesting direction for
research.

## References

  * Ailon et al. [2021] Nir Ailon, Omer Leibovitch, and Vineet Nair.  Sparse linear networks with a fixed butterfly structure: theory and practice.  In _Uncertainty in Artificial Intelligence_ , pages 1174–1184. PMLR, 2021.
  * Bartlett et al. [2022] Peter Bartlett, Piotr Indyk, and Tal Wagner.  Generalization bounds for data-driven numerical linear algebra.  In _Conference on Learning Theory_ , pages 2013–2040. PMLR, 2022\.
  * Chierichetti et al. [2017] Flavio Chierichetti, Sreenivas Gollapudi, Ravi Kumar, Silvio Lattanzi, Rina Panigrahy, and David P Woodruff.  Algorithms for $\ell_{p}$ low-rank approximation.  In _International Conference on Machine Learning_ , pages 806–814. PMLR, 2017.
  * Clarkson and Woodruff [2009] Kenneth L Clarkson and David P Woodruff.  Numerical linear algebra in the streaming model.  In _Proceedings of the forty-first annual ACM symposium on Theory of computing_ , pages 205–214, 2009.
  * Clarkson and Woodruff [2017] Kenneth L Clarkson and David P Woodruff.  Low-rank approximation and regression in input sparsity time.  _Journal of the ACM (JACM)_ , 63(6):1–45, 2017\.
  * Cohen et al. [2015] Michael B Cohen, Sam Elder, Cameron Musco, Christopher Musco, and Madalina Persu.  Dimensionality reduction for k-means clustering and low rank approximation.  In _Proceedings of the forty-seventh annual ACM symposium on Theory of computing_ , pages 163–172, 2015.
  * Gupta and Roughgarden [2016] Rishi Gupta and Tim Roughgarden.  A pac approach to application-specific algorithm selection.  In _Proceedings of the 2016 ACM Conference on Innovations in Theoretical Computer Science_ , pages 123–134, 2016.
  * Gupta and Roughgarden [2020] Rishi Gupta and Tim Roughgarden.  Data-driven algorithm design.  _Communications of the ACM_ , 63(6):87–94, 2020\.
  * Halko et al. [2011] Nathan Halko, Per-Gunnar Martinsson, and Joel A Tropp.  Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions.  _SIAM review_ , 53(2):217–288, 2011.
  * Imamoglu et al. [2018] Nevrez Imamoglu, Yu Oishi, Xiaoqiang Zhang, Guanqun Ding, Yuming Fang, Toru Kouyama, and Ryosuke Nakamura.  Hyperspectral image dataset for benchmarking on salient object detection.  In _2018 Tenth international conference on quality of multimedia experience (qoMEX)_ , pages 1–3. IEEE, 2018.
  * Indyk [2006] Piotr Indyk.  Stable distributions, pseudorandom generators, embeddings, and data stream computation.  _Journal of the ACM (JACM)_ , 53(3):307–323, 2006\.
  * Indyk et al. [2019] Piotr Indyk, Ali Vakilian, and Yang Yuan.  Learning-based low-rank approximations.  _Advances in Neural Information Processing Systems_ , 32, 2019.
  * Indyk et al. [2021] Piotr Indyk, Tal Wagner, and David Woodruff.  Few-shot data-driven algorithms for low rank approximation.  _Advances in Neural Information Processing Systems_ , 34:10678–10690, 2021.
  * Koch and Lubich [2007] Othmar Koch and Christian Lubich.  Dynamical low-rank approximation.  _SIAM Journal on Matrix Analysis and Applications_ , 29(2):434–454, 2007.
  * Liberty et al. [2007] Edo Liberty, Franco Woolfe, Per-Gunnar Martinsson, Vladimir Rokhlin, and Mark Tygert.  Randomized algorithms for the low-rank approximation of matrices.  _Proceedings of the National Academy of Sciences_ , 104(51):20167–20172, 2007.
  * Liu et al. [2022] Cuiyu Liu, Chuanfu Xiao, Mingshuo Ding, and Chao Yang.  Tensor-based sketching method for the low-rank approximation of data streams.  _arXiv preprint arXiv:2209.14637_ , 2022.
  * Liu et al. [2020] Simin Liu, Tianrui Liu, Ali Vakilian, Yulin Wan, and David P Woodruff.  Learning the positions in countsketch.  _arXiv preprint arXiv:2007.09890_ , 2020.
  * Sakaue and Oki [2022] Shinsaku Sakaue and Taihei Oki.  Improved generalization bound and learning of sparsity patterns for data-driven low-rank approximation.  _arXiv preprint arXiv:2209.08281_ , 2022.
  * Sarlos [2006] Tamas Sarlos.  Improved approximation algorithms for large matrices via random projections.  In _2006 47th annual IEEE symposium on foundations of computer science (FOCS’06)_ , pages 143–152. IEEE, 2006.
  * Woolfe et al. [2008] Franco Woolfe, Edo Liberty, Vladimir Rokhlin, and Mark Tygert.  A fast randomized algorithm for the approximation of matrices.  _Applied and Computational Harmonic Analysis_ , 25(3):335–366, 2008.