Quantum-inspired algorithms in practice

I Introduction

A driving force for studying quantum computing is the conviction that quantum algorithms can solve some problems more efficiently than classical methods. But the boundary between classical and quantum computing changes constantly: new classical and quantum algorithms disrupt the landscape, continuously updating the border that separates these two paradigms. This interplay is fruitful; it allows classical and quantum algorithm developers to feed off each other’s innovations and advance our understanding of the limits of computation. Indeed, several classical algorithms have been reported that rely at least partially on developments in quantum computing Narayanan and Moore (1996); Han and Kim (2002); da Cruz et al. (2007); Valiant (2008); Rønnow et al. (2014); Barak et al. (2015); Chakhmakhchyan et al. (2017); Tang (2018a); Chia et al. (2019).

Quantum algorithms for linear algebra are a flagship application of quantum computing, particularly due to their relevance in machine learning

Harrow et al. (2009); Lloyd et al. (2014); Rebentrost et al. (2014); Kerenidis and Prakash (2016); Biamonte et al. (2017). These algorithms typically scale polylogarithmically with dimension, which, at the time they were reported, implied an asymptotic exponential speedup compared to state-of-the-art classical methods. For this reason, significant interest has been generated in the dequantization approach that led to breakthrough quantum-inspired classical algorithms for linear algebra problems with sublinear complexity Tang (2018a, b); Gilyén et al. (2018); Chia et al. (2018)

. These dequantized algorithms work for general low-rank matrices, whereas quantum computers still exhibit an exponential speedup over all known classical algorithms for sparse, full-rank matrix problems, including the quantum Fourier transform, eigenvector and eigenvalue analysis, linear systems, and others. Dequantized algorithms for such problems – high-rank matrix inversion, for example – would imply that classical computers can efficiently simulate quantum computers, i.e., BQP = P, which is not currently considered to be likely.

The proof techniques used to determine the complexity of quantum-inspired algorithms lead to bounds that suggest that runtimes might be prohibitively large for practical applications: the proven complexity of the linear systems algorithm is $~ O (κ^{16} k^{6} ∥ A ∥_{F}^{6} / ε^{6})$ Gilyén et al. (2018), while for recommendation systems it is $~ O (k^{12} / ε^{12})$ Tang (2018a)¹¹1We have taken $σ = O (∥ A ∥_{F} / \sqrt{k})$ as stated in Tang (2018a)

and ignored scaling with respect to the failure probability

η

.. Here,

A

is the input matrix,

∥ A ∥_{F}

is the Frobenius norm of

A

κ

is the condition number,

k

is the rank, and

ε

is the precision of the algorithm. These results raise a few immediate questions: are these complexity bounds the consequence of proof techniques, or do they reflect a fundamental limitation on the capabilities of these algorithms? How do these algorithms actually perform in practice?

We address these questions by analyzing, implementing, and testing quantum-inspired algorithms for linear algebra. First, we perform a theoretical analysis aimed at identifying potential bottlenecks in their practical implementation. This allows us to anticipate the regimes where quantum-inspired algorithms may outperform previously-known classical methods. We then test the algorithms on artifical examples based on random matrices, where it is possible to control the dimension, rank, and condition number of the input matrix. We conclude by testing the algorithms on problems of practical interest. Specifically, we run the quantum-inspired algorithm for linear systems of equations applied to portfolio optimization on stocks from the S&P 500, and analyze the algorithm for recommendation systems on the MovieLens dataset Harper and Konstan (2016). Based on our analysis and tests, we find that, provided that the input matrices have small rank and condition number, quantum-inspired algorithms can perform well in practice: they provide good approximations in reasonable times, even for very large-dimensional problems. This shows that quantum-inspired algorithms can have significantly smaller runtimes than their worst-case complexity bounds would suggest, although we find evidence that the dependence on the error $ε$ may be tight for the linear systems algorithm. For problems whose input matrices have larger ranks and condition numbers, the algorithms struggle noticeably to produce good results.

In the following sections, we begin by providing a brief description of quantum-inspired algorithms for linear algebra. This high-level summary may prove useful for those wishing to become more familiar with these techniques. We then analyze each of the main steps of the algorithm in further detail, focusing on identifying potential practical bottlenecks. We continue by implementing the algorithms and benchmarking their performance on artificially-generated problems and on real-world data. We conclude by summarizing the implications of our analysis.

Ii Quantum-inspired algorithms for linear algebra

In this section, we give an overview of quantum-inspired algorithms for linear systems of equations and recommendation systems. More details can be found in Refs. Tang (2018a, b); Gilyén et al. (2018); Chia et al. (2018)

. A unique feature of our description is that we view both algorithms as specific examples of a more general method to sample from vectors expressed in terms of the singular value decomposition (SVD) of an input matrix. Henceforth, capital letters are employed to denote matrices and bold letters to denote column vectors. For example, a linear system of equations is written as

A x = b

, where

A \in R^{m \times n}

x = (x_{1}, x_{2}, \dots, x_{n})^{T}

, and

b = (b_{1}, b_{2}, \dots, b_{m})^{T}

We consider problems of the following form: Given an $m \times n$ matrix $A \in R^{m \times n}$ with SVD

A = k \sum ℓ = 1 σ_{ℓ} u^{(ℓ)} {v^{(ℓ)}}^{T},

(1)

the goal is to sample entries of the $n$ -dimensional vector

x = k \sum ℓ = 1 λ_{ℓ} v^{(ℓ)},

(2)

with respect to the length-squareprobability distribution $p_{x} (i) = x_{i}^{2} / ∥ x ∥^{2}$ . The coefficients $λ_{ℓ}$ depend on the matrix $A$ and possibly on other inputs. For linear systems of equations, $x$ is the solution vector given by $x = A^{+} b$ , with $A^{+}$ the Moore-Penrose pseudoinverse. As discussed in Ref. Gilyén et al. (2018), the coefficients $λ_{ℓ}$ of Eq. (2) are given by the inner product

λ_{ℓ} = \frac{1}{σ_{ℓ}^{2}} ⟨ v^{(ℓ)}, A^{T} b ⟩ .

(3)

Similarly, for recommendation systems, $A$ is the preference matrix, whose entries $A_{i j}$ denote the rating that user $i$ has given to product $j$ . In this case, the solution vector $x$ is the $i$ -th row of a low-rank approximation of $A$ , which contains all the preferences of user $i$ . The coefficients are given by

λ_{ℓ} = ⟨ A_{i}^{T}, v^{(ℓ)} ⟩,

(4)

where $A_{i}$ is the $i$ -th row of $A$ Tang (2018a).

Quantum-inspired algorithms assume that the entries of $A$ are given in a way that allows length-square sampling to be performed. This can be accomplished, for example, if the entries of $A$ are stored in a data structure of the form proposed by Kerenidis and Prakash (2016), which stores both entries of vectors and their sub-norms. Length-square sampling allows us to preferentially sample the rows of $A$ with large norm, and the large entries within each row. Given access to length-squared sampling, quantum-inspired algorithms solve low-rank linear algebra problems in three main steps:

Approximate SVD: The Frieze-Kannan-Vempala (FKV) algorithm Frieze et al. (2004) is used to compute approximate singular values ${~ σ}_{ℓ}$ and approximate right singular vectors ${~ v}^{(ℓ)}$ . The singular vectors are not calculated in their entirety; instead, a description is found allowing efficient computation of any given entry $v_{i}^{(ℓ)}$ .
Coefficient estimation: Based on the results of step 1, the coefficients ${~ λ}_{ℓ}$
are approximated using Monte Carlo estimation techniques.
Sampling solution vectors: Using Eq. (2) and the results from steps 1 and 2, rejection sampling is employed to perform length-square sampling from the approximate vectors $~ x = \sum_{ℓ = 1}^{k} {~ λ}_{ℓ} {~ v}^{(ℓ)}$ .

Before describing each of these steps in more detail, it is important to understand that the innovations of quantum-inspired algorithms are steps 2 and 3. In these steps, coefficient estimation and sampling from the solution vectors is performed in time $O (poly (k, κ, ε, log n, log m))$ . By contrast, the strategy originally suggested in Ref. Frieze et al. (2004) is simply to compute coefficients and solution vectors directly from the approximate SVD of step 1, which requires time $O (k n)$ . Asymptotically, quantum-inspired algorithms thus achieve an exponential asymptotic speedup in dimension, from $O (n, m)$ to $polylog (n, m)$ , at the cost of a polynomial dependency on other parameters. In practice, the direct calculation can be done extremely fast since it scales only linearly with dimension, meaning that quantum-inspired algorithms require very large-dimensional problems before it becomes preferable to employ their sampling techniques.

ii.1 Approximate SVD

The Frieze-Kannan-Vempala (FKV) algorithm Frieze et al. (2004) is a method to compute low-rank approximations of matrices. The FKV algorithm is just a particular example among several randomized methods for computing approximate matrix decompositions Drineas et al. (2006); Sarlos (2006); Woolfe et al. (2008); Rokhlin and Tygert (2008); Rokhlin et al. (2009); Martinsson et al. (2011), see Ref. Halko et al. (2011) for a survey of these techniques. It is an interesting question whether other randomized methods for approximate matrix decompositions can be used to improve quantum-inspired algorithms. For concreteness, in this work we focus on analyzing the FKV algorithm.

The main idea behind FKV is the following: instead of performing a singular value decomposition of the large input matrix $A \in R^{m \times n}$ , a smaller matrix $C \in R^{r \times c}$ is constructed by sampling $r$ rows and $c$ columns of $A$ . The only restriction on $r$ and $c$ is that they have to be larger than the rank of $A$ : in particular, they can be independent of $m$ and $n$ , i.e., not even $O (log m)$ , $O (log n)$ . The matrix $C$ captures the important features of $A$ , making it possible to perform an SVD of the smaller matrix $C$ to retrieve information about the singular values and vectors of $A$ . Since the rank of $C$ is bounded by its dimension, this approach is only possible when $A$ is low-rank or when a low-rank approximation of $A$ is sufficient.

Figure 1: Schematic representation of the FKV algorithm. The input is an $m \times n$ matrix $A$ , which in this example has ten rows and columns. After performing one pass through the matrix, it is possible to perform length-square sampling of rows and columns. A total of $r$ rows are sampled from this length-square distribution, then renormalized and stacked together to construct a new $r \times n$ matrix $R$ . In this example, we set $r = 4$ and obtain the row indices $(i_{1}, i_{2}, i_{3}, i_{4}) = (2, 8, 3, 6)$ . Once the matrix $R$ has been built, $c$ column indices are drawn by repeatedly choosing rows of $R$ uniformly at random, then sampling from the corresponding length-square distribution of column indices. In this example, we set $c = 4$ and obtain the column indices $(j_{1}, j_{2}, j_{3}, j_{4}) = (3, 8, 1, 3)$ . Note that repeated indices are allowed in the algorithm. The final matrix $C$ is built by renormalizing these columns of $R$ and grouping them together. The singular values of $C$ approximate those of $A$ , while the singular vectors of $C$ can be used to reconstruct approximations to the singular vectors of $A$ .

In more detail, the first step of FKV is to perform a pass through the matrix $A$ and compute the Frobenius norm of each row $∥ A_{i} ∥_{F}$ , as well as the norm of the matrix $∥ A ∥_{F}$ . These quantitites define a length-square distribution over rows: $p (i) = ∥ A_{i} ∥_{F}^{2} / ∥ A ∥_{F}^{2}$ . For each row $i$ , a length-square distribution over its entries is also computed, $q_{i} (j) = ∥ A_{i j} ∥_{F}^{2} / ∥ A_{i} ∥_{F}^{2}$ . To sample from these distributions, it is possible to employ a specialized data structure as suggested in Refs. Frieze et al. (2004); Kerenidis and Prakash (2016); Tang (2018a) that allows sampling in time that is logarithmic in dimension. The algorithm proceeds as follows:

Sample $r$ row indices $i_{1}, i_{2}, \dots, i_{r}$ from the row distribution $p (i)$ . For each row index, select the row $A_{i_{s}}$ and renormalize it as $R_{i_{s}} = \frac{∥ A ∥_{F}}{\sqrt{r} ∥ A_{i_{s}} ∥} A_{i_{s}}$ . This defines an $r \times n$ matrix $R$ consisting of the rescaled rows $R_{i_{s}}$ .
Select an index $s \in {1, 2, \dots, r}$ uniformly at random, then sample a column index $j$ from the column distribution $q_{i_{s}} (j)$ . Repeat this a total of $c$ times to obtain column indices $j_{1}, j_{2}, \dots, j_{c}$ . For each column index, select the column $R_{\cdot, j_{t}}$ and renormalize it as $C_{\cdot, j_{t}} = \frac{∥ A ∥_{F}}{\sqrt{c} ∥ R_{\cdot, j_{t}} ∥} R_{\cdot, j_{t}}$ . This defines a matrix $C$ consisting of the rescaled columns $C_{\cdot, j_{t}}$ .
Compute the singular value decomposition of $C$ .

These steps are illustrated in Fig. 1. The FKV algorithm performs a useful tradeoff: it is possible to perform SVDs of a significantly smaller matrix $C$ at the cost of obtaining only approximations of the singular values and vectors of the original matrix $A$ . We denote by ${~ σ}_{ℓ}$ the singular values of $C$ and by $ω^{(ℓ)}$ its left singular vectors. The FKV algorithm directly provides approximations of singular values, ${~ σ}_{ℓ} \approx σ_{ℓ}$ for $ℓ = 1, 2, \dots, k$ . Approximate right singular vectors ${~ v}^{(ℓ)}$ of $A$ are obtained using the expression

{~ v}^{(ℓ)} = \frac{1}{{~ σ}_{ℓ}} R^{T} ω^{(ℓ)},

(5)

while approximate left singular vectors ${~ u}^{(ℓ)}$ are given by

{~ u}^{(ℓ)} = A (\frac{1}{{~ σ}_{ℓ}^{2}} R^{T} ω^{(ℓ)}) .

(6)

In existing classical methods, these vectors are computed directly, whereas in quantum-inspired algorithms it suffices to query their entries.

The quality of the FKV approximations depends on the number of sampled rows $r$ and columns $c$ . The challenge is to find values of $r$ and $c$ such that sufficiently good approximations are obtained, while also ensuring that the complexity of computing the SVD of matrix $C$ is significantly less than for the original input matrix $A$ , i.e., such that $r ≪ m$ and $c ≪ n$ . Since $r, c$ can be much smaller than $m, n$ , a substantial runtime reduction can potentially be achieved.

FKV is the core of quantum-inspired algorithms for linear algebra. Without the use of randomized methods, computing the SVD of an $m \times n$ matrix $A$ requires $O (min {m^{2} n, m n^{2}})$ time using naive matrix multiplication. From that point onwards, computing coefficients $λ$ and the solution vector $x$ takes only linear time in $m, n$ . In FKV, we instead calculate the SVD of an $r \times c$ matrix, which requires $O (min {r^{2} c, r c^{2}})$ time. When $r, c$ are much smaller than $m, n$ , this where the real savings of quantum-inspired algorithms happen. As explained previously, the subsequent steps of coefficient estimation and sampling from the solution vector are used to reduce a linear runtime to polylogarithmic runtime. As we discuss later, for practical problem sizes that are far from the asymptotic limit, it is preferable to compute coefficients and solution vectors explicitly by employing Eqs. (5) and (6) starting from the approximate SVD of matrix $C$ in the FKV algorithm.

ii.2 Coefficient estimation

The coefficients $λ$ appearing in Eq. (2) are inner products between vectors, multiplied by a power of the singular values in the case of linear systems. Quantum-inspired algorithms compute approximations of these coefficients using Monte Carlo estimation. In general, a coefficient $λ$ is given by

λ = ⟨ y, z ⟩ = n \sum i = 1 y_{i} z_{i},

(7)

for some appropriate vectors $y, z$ , as in Eqs. (3) and (4). The strategy to estimate this inner product is to perform length-square sampling from one of the vectors. Without loss of generality, we take this vector to be $y$

. Define the random variable

χ

that takes values

χ_{i} = y_{i} z_{i} / p_{x} (i)

, where the indices

i

are sampled from the length-square distribution

p_{y} (i) = y_{i}^{2} / ∥ y ∥^{2}

. The expectation value of the random variable satisfies

E (χ) = n \sum i = 1 \frac{y_{i} z_{i}}{p_{y} (i)} p_{y} (i) = ⟨ y, z ⟩ = λ .

(8)

Similarly, the second moment is

E (χ^{2}) = n \sum i = 1 {(\frac{y_{i} z_{i}}{p_{y} (i)})}^{2} p_{y} (i) = ∥ y ∥^{2} ∥ z ∥^{2},

(9)

and $σ_{χ}^{2} = E (χ^{2}) - E (χ)^{2}$

is the variance of

χ

. The strategy to estimate coefficients

λ

is to draw

N

samples

χ^{(1)}, χ^{(2)}, \dots, χ^{(N)}

from

χ

and compute the unbiased estimator

^λ = \frac{1}{N} \sum_{j = 1}^{N} χ^{(j)} \approx λ

. This constitutes a form of importance sampling in Monte Carlo estimation: large entries of

y

are preferably selected since they contribute more significantly to the inner product. The error in the estimation is quantified by the ratio between the standard deviation and the mean of the estimator:

ϵ := \sqrt{Var (^λ)} / | λ |

. The variance of the estimator is

Var (^λ) = σ_{χ}^{2} / N

, leading to a precision

ϵ = σ_{χ} / (| λ | \sqrt{N})

in the estimator. This implies that the number of samples

N

needed to, with high probability, achieve a precision

ϵ

	$N =$	$\frac{σ_{χ}^{2}}{λ^{2} ϵ^{2}} = \frac{1}{ϵ^{2}} [\frac{∥ y ∥^{2} ∥ z ∥^{2}}{⟨ y, z ⟩^{2}} - 1]$		(10)
	$=$	$O (\frac{∥ y ∥^{2} ∥ z ∥^{2}}{ϵ^{2} ⟨ y, z ⟩^{2}}) = O (\frac{1}{ϵ^{2} {cos}^{2} θ}),$		(11)

where $θ$ is the angle between the vectors $y$ and $z$ .

ii.3 Complexity of coefficient estimation

Quantum-inspired algorithms differ from the FKV algorithm because they estimate the coefficients $λ$ rather than calculating them directly from the approximate SVD of the input matrix. Therefore, it is worthwhile to examine the complexity of coefficient estimation step to determine the regime where quantum-inspired algorithms become beneficial. For $n$ -dimensional random vectors, the expectation value of the angle between two vectors satisfies $E [\frac{1}{{cos}^{2} θ}] = n$ , so we can anticipate that the number of samples needed may scale linearly with dimension. If $\frac{1}{{cos}^{2} θ} = poly (m, n)$ , the algorithm no longer has the desired polylogarithmic complexity. In fact, as we show in Appendix A, regardless of what coefficient estimation strategy is employed, in a worst-case setting, inner products of vectors cannot be approximated in sublinear time. However, for matrices of low rank and condition number, as we show below, it is possible to achieve good estimation of coefficients.

For linear systems of equations, the coefficients $λ_{ℓ}$ are given by

λ_{ℓ} = \frac{1}{σ_{ℓ}^{2}} ⟨ v^{(ℓ)}, A^{T} b ⟩ = \frac{1}{σ_{ℓ}} ⟨ u^{(ℓ)}, b ⟩ .

(12)

In the quantum-inspired algorithm, we do not have direct access to the length-square distributions of either $A^{T} b$ or $v^{(ℓ)}$ . Instead, the strategy in Ref. Gilyén et al. (2018) is to express the coefficients as $λ_{ℓ} = \frac{1}{σ_{ℓ}^{2}} Tr (A^{T} b {v^{(ℓ)}}^{T})$ and sample the entries $A_{i j}$ of $A$ . In this case the the number of samples needed in the estimation is Gilyén et al. (2018):

N

= \frac{1}{ϵ^{2}} (\frac{∥ v^{(ℓ)} ∥^{2} ∥ A ∥_{F}^{2} ∥ b ∥^{2}}{⟨ v^{(ℓ)}, A^{T} b ⟩^{2}} - 1) = O (\frac{1}{ϵ^{2}} \frac{∥ σ ∥^{2}}{σ_{ℓ}^{2}} \frac{∥ b ∥^{2}}{β_{ℓ}^{2}}),

(13)

where we have defined $σ := (σ_{1}, σ_{2}, \dots, σ_{k})$ , $β_{ℓ} := ⟨ u^{(ℓ)}, b ⟩$ such that $b = \sum_{ℓ = 1}^{k} β_{ℓ} u^{(ℓ)}$ , and we have used Eq. (12) as well as the identity $∥ A ∥_{F}^{2} = ∥ σ ∥^{2}$ . In the special cases of recommendation systems and portfolio optimization, the coefficients are proportional to an inner product $⟨ A_{i}^{T}, v^{(ℓ)} ⟩$ . Writing $A_{i} = \sum_{ℓ^{'} = 1}^{k} σ_{ℓ^{'}} u_{i}^{(ℓ^{'})} {v^{(ℓ^{'})}}^{T}$ we have

⟨ A_{i}^{T}, v^{(ℓ)} ⟩ = k \sum ℓ^{'} = 1 σ_{ℓ^{'}} u_{i}^{(ℓ^{'})} ⟨ v^{(ℓ^{'})}, v^{(ℓ^{'})} ⟩ = σ_{ℓ} u_{i}^{(ℓ)} .

(14)

The number of samples then satisfies

N

= \frac{1}{ϵ^{2}} (\frac{∥ A_{i}^{T} ∥^{2} ∥ v^{(ℓ)} ∥^{2}}{⟨ A_{i}^{T}, v^{(ℓ)} ⟩^{2}} - 1) = O ⎛ ⎜ ⎝ \frac{\sum_{ℓ^{'} = 1}^{k} (σ_{ℓ^{'}} u_{i}^{(ℓ^{'})})^{2}}{ϵ^{2} σ_{ℓ}^{2} {u_{i}^{(ℓ)}}^{2}} ⎞ ⎟ ⎠ .

(15)

Define the $k$ -dimensional vectors $σ := (σ_{1}, σ_{2}, \dots, σ_{k})^{T}$ , $μ := (u_{i}^{(1)}, u_{i}^{(2)}, \dots, u_{i}^{(k)})^{T}$ and their Hadamard product $ν = σ \circ μ$ with entries $ν_{ℓ} = σ_{ℓ} u_{i}^{(ℓ)}$ . It follows that

N = O (\frac{1}{ϵ^{2}} \frac{∥ ν ∥^{2}}{ν_{ℓ}^{2}}) .

(16)

In all cases, the number of samples depends on the ratio between the norm squared of a $k$ -dimensional vector and the square of the $ℓ$ -th entry of that vector. We now determine the worst-case scaling for these ratios. We consider the ratio $∥ σ ∥^{2} / σ_{ℓ}^{2}$ for concreteness; the same will apply for $b$ and $ν$ . The largest ratio occurs for the smallest singular value $σ_{min}$ and we have

\frac{∥ σ ∥^{2}}{σ_{min}^{2}} = \frac{∥ σ ∥^{2}}{} σ_{max}^{2} \frac{σ_{max}^{2}}{σ_{min}^{2}} = \frac{∥ σ ∥^{2}}{σ_{max}^{2}} κ^{2} .

(17)

Similarly, this ratio is largest when $σ_{max}$ is as small as possible. This occurs when the largest $k - 1$ singular values are equal to each other. In that case, defining $σ_{ℓ}^{'} := σ_{ℓ} / ∥ σ ∥$ such that $∥ σ^{'} ∥ = 1$ we have

	$∥ σ^{'} ∥^{2} = (σ_{max}^{'})^{2} (k - 1 + \frac{1}{κ}) = 1$		(18)
	$\Rightarrow (σ_{max}^{'})^{2} = \frac{1}{k - 1 + \frac{1}{κ}} = O (1 / k),$		(19)

and therefore

\frac{∥ σ ∥^{2}}{σ_{max}^{2}} κ^{2} = \frac{1}{(σ_{max}^{'})^{2}} κ^{2} = O (k κ^{2}) .

(20)

Defining $κ_{β} := β_{max} / β_{min}$ and $κ_{ν} := ν_{max} / ν_{min}$ , where $β_{max/min}$ and $ν_{max/min}$ are the largest and smallest entries of their corresponding vectors, we conclude that the number of samples required to estimate coefficients in quantum-inspired algorithms for linear systems scales as

N = O (\frac{k^{2} κ^{2} κ_{β}^{2}}{ϵ^{2}}),

(21)

and for recommendation systems as

N = O (\frac{k κ_{ν}^{2}}{ϵ^{2}}) .

(22)

These asymptotic formulas – or more precisely Eqs. (13) and (15) – can be used to estimate the problem sizes for which the complexity of coefficient estimation is smaller than a direct calculation. As an illustration, in linear systems, setting $k = κ = κ_{β} = 100$ means that an order of $N = 10^{16}$ samples are needed to estimate the coefficients with an error of $ϵ = 10^{- 2}$ . Finally, it is important to understand that a large number of samples can improve the precision of the estimation, but not its accuracy: if the singular values and vectors obtained from FKV have large errors, the expectation value of the estimator will not coincide with the actual value of the coefficient.

ii.4 Sampling solution vectors

The approximate SVD and coefficient estimation steps respectively provide approximate singular values ${~ σ}_{ℓ}$ and approximate coefficients ${~ λ}_{ℓ}$ . The approximate right singular vectors of $A$ are not constructed explicitly, but are implicitly defined by Eq. (5), namely ${~ v}^{(ℓ)} = \frac{1}{{~ σ}_{ℓ}} R^{T} w^{(ℓ)}$ . This allows an implicit calculation of the approximate solution vector $~ x = R^{T} w$ , where

w := k \sum ℓ = 1 \frac{{~ λ}_{ℓ}}{{~ σ}_{ℓ}} w^{(ℓ)} .

(23)

The challenge is to sample from these vectors using only query access to the entries of $w$ and $R$ . Quantum-inspired algorithms achieve this using a well-known technique known as rejection sampling. The steps are as follows Tang (2018a):

Sample a row index $i$ uniformly at random.
Sample a column index $j$ from the length-square distribution $q_{i} (j) = | R_{i j} |^{2} / ∥ R_{i} ∥^{2}$ .
Output $j$ with probability $\frac{| ⟨ w, R_{\cdot, j} ⟩ |^{2}}{∥ R_{\cdot, j} ∥^{2} ∥ w ∥^{2}}$ , or sample $j$ again.

The expected number of repetitions before outputting a column index $j$ is $r ∥ w ∥^{2} / (\sum_{j = 1}^{n} | ⟨ w, R_{\cdot, j} ⟩ |^{2}) = r ∥ w ∥^{2} / ∥ ~ x ∥^{2}$ . The steps of rejection sampling can be performed quickly. The overhead arises from calculating the inner product $⟨ w, R_{\cdot, j} ⟩$ , which requires $O (r k)$ time; significantly less than computing the approximate SVD. As we illustrate in the following sections, the number of expected repetitions is usually moderate in practice, so sampling from the approximate solution vector can be performed in comparatively less time than other steps in the algorithm.

Iii Numerical benchmarking

In the previous section we performed a theoretical analysis of quantum-inspired algorithms. Here, we complement that analysis by benchmarking the algorithms on specific problems. We examine the performance of each step of the algorithm individually, revealing the ones that contribute most to the overall runtime and which ones lead to the most significant source of errors. As an initial benchmark, we study the algorithm for linear systems of equations on artificial problems of extremely large dimension. We design matrices such that performing length-square sampling from its rows and columns, as well as all other steps of the quantum-inspired algorithm can be done using only query access to its entries. This allows to study the performance of the algorithms in their intended asymptotic regime, without handling large datasets directly. We then test the algorithm for linear systems on Gaussian random matrices. This choice is made to allow full control over properties of the input matrices: dimension, rank, and condition number, while collecting large amounts of statistics. Afterwards, we apply the algorithms to portfolio optimization and movie recommendations. All algorithms were implemented in Python and the source code is available at https://github.com/XanaduAI/quantum-inspired-algorithms.

iii.1 High-dimensional problems

Quantum-inspired algorithms for linear algebra are aimed at tackling extremely large datasets, but it is very challenging to test them in such scenarios. Instead, we design a class of problems where it is possible to run quantum-inspired algorithms using only black-box access to the entries of the input matrix and vector. We focus on the algorithm for linear systems and report the main formulas, whose derivations can be found in Appendix B.

Define the matrices

V^{(ℓ)} = \frac{1}{2^{n}} \sum y, z \in {0, 1}^{n} (- 1)^{x^{(ℓ)} \cdot (y \oplus z)} e^{(y)} {e^{(z)}}^{T},

(24)

where $x, y, z$ are $n$ -bit strings, the vectors $e^{(y)}$ form a basis of the $2^{n}$ -dimensional vector space, and the $x^{(ℓ)}$ are distinct $n$ -bit strings for $ℓ = 1, 2, \dots, k$ . We define the input matrix as

A :=

k \sum ℓ = 1 σ_{ℓ} V^{(ℓ)} = \frac{1}{2^{n}} \sum y, z \in {0, 1}^{n} a_{y, z} e^{(y)} {e^{(z)}}^{T},

(25)

where $σ_{1}, σ_{2}, \dots, σ_{k}$ are the singular values and we have implicitly defined $a_{y, z} = \sum_{ℓ = 1}^{k} σ_{ℓ} (- 1)^{x^{(ℓ)} \cdot (y \oplus z)}$ . The right and left singular vectors of $A$ are the same, given by $v^{(ℓ)} = \frac{1}{\sqrt{2^{n}}} \sum_{y \in {0, 1}^{n}} (- 1)^{x^{(ℓ)} \cdot y} e^{(y)}$ . As we show in the Appendix B, by construction, all rows of $A$ have the same norm, therefore length-square sampling can be done by choosing rows uniformly at random. For length-square sampling of columns, we use rejection sampling: given a row $y$ , we select a column index $z$ uniformly at random and accept it with probability

p (z) = \frac{| a_{y, z} |^{2}}{{max}_{z^{'}} | a_{y, z^{'}} |^{2}} .

(26)

The optimization is trivial since ${max}_{z^{'}} | a_{y, z^{'}} |^{2} = \sum_{ℓ = 1}^{k} (σ_{ℓ})^{2}$ . Once row indices $i_{1}, \dots, i_{r}$ and column indices $j_{1}, \dots, j_{c}$ have been sampled, the entries of the $r \times c$ matrix $C$ in FKV can be expressed as

C_{s, t} = \frac{1}{\sqrt{2^{k} c}} \frac{∥ a ∥}{\sqrt{\sum_{s^{'} = 1}^{r} (a_{i_{s^{'}}, j_{t}})^{2}}} a_{i_{s}, j_{t}},

(27)

where $a = (a_{1}, a_{2}, \dots, a_{2^{k}})^{T}$ is a vector containing all the possible values of the entries $a_{y, z}$ . The approximate singular values and vectors of $A$ can be reproduced from the SVD of $C$ as detailed in previous sections. To estimate the coefficient $λ_{ℓ} = \frac{1}{σ_{ℓ}^{2}} ⟨ v^{(ℓ)}, A^{T} b ⟩$ , we sample the random variable

χ_{y, z} = \frac{∥ a ∥^{2}}{2^{k}} \frac{(- 1)^{x^{(ℓ)} \cdot z}}{\sum_{ℓ = 1}^{k} σ_{ℓ} (- 1)^{x^{(ℓ)} \cdot (y \oplus z)}} k \sum ℓ = 1 β_{ℓ} (- 1)^{x^{(ℓ)} \cdot y},

(28)

and use its sample mean as an estimator for the coefficients. Finally, the approximate solution vector $~ x$ is given by

~ x = \frac{1}{\sqrt{2^{n}}} k \sum ℓ = 1 {~ λ}_{ℓ} {~ v}^{(ℓ)} .

(29)

Any entry of the approximate solution vector can be computed efficiently from the entries of the approximate right singular vectors ${~ v}^{(ℓ)}$ , which can be obtained from the approximate SVD of $C$ .

The crucial property of this construction is that the only dependency on dimension comes from the sampling of rows and columns, which can be done efficiently even for extremely high-dimensional problems. The limitations arise in computing the matrix $C$ , which takes $O (r^{2} c)$ time, and from computing the $2^{k}$ different matrix elements. This limits the values of $r, c$ and $k$ that can be explored in this approach. We fix the input matrices to have an extremely large dimension of $2^{50} \approx 10^{15}$ . We set $r = c = 150$ , the largest that could be handled in a reasonable amount of time. This means that we are aiming to approximate the SVD of the input matrix by using a smaller matrix whose dimension differs by thirteen orders of magnitude. The bit strings $x^{(ℓ)}$ that define each singular vector are chosen uniformly at random. We study three examples with $k = κ = κ_{β} = 3, 5, 10$ , employing $N = 10^{4}$ samples in the estimation of the coefficients. We report the average and standard deviation of the errors $η_{σ} = \frac{1}{k} \sum_{ℓ = 1}^{k} \frac{| {~ σ}_{ℓ} - σ_{ℓ} |}{σ_{ℓ}}$ and $η_{λ} = \frac{1}{k} \sum_{ℓ = 1}^{k} \frac{| {~ λ}_{ℓ} - λ_{ℓ} |}{| λ_{ℓ} |}$ over ten repetitions of the algorithms. For vectors, we evaluate the first $L$ entries and similarly report the average and standard deviation of the errors $η_{v} = \frac{1}{k L} \sum_{i = 1}^{L} \sum_{ℓ = 1}^{k} \frac{| v_{i}^{(ℓ)} - {~ v}_{i}^{(ℓ)} |}{| v_{i}^{(ℓ)} |}$ and $η_{x} = \frac{1}{k L} \sum_{i = 1}^{L} \sum_{ℓ = 1}^{k} \frac{| x_{i} - {~ x}_{i} |}{| x_{i} |}$ . The results are summarized in Table 1.

	Errors				•
Case study	$η_{σ}$	$η_{v}$	$η_{λ}$	$η_{x}$	Runtime (hours)
$k = κ = κ_{β} = 3$	$0.011 \pm 0.012$	$0.124 \pm 0.071$	$0.285 \pm 0.276$	$0.414 \pm 0.154$	4.3
$k = κ = κ_{β} = 5$	$0.129 \pm 0.016$	$0.212 \pm 0.070$	$0.530 \pm 0.616$	$1.235 \pm 0.615$	8.6
$k = κ = κ_{β} = 10$	$0.626 \pm 0.018$	$1.619 \pm 0.264$	$1.193 \pm 1.720$	$4.138 \pm 2.321$	29.0

Table 1: Errors in estimating singular values, singular vectors, coefficients, and solution vectors for input matrices of dimension

n, m = 2^{50}

. The coefficients were estimated using

N = 10^{4}

samples and the errors in the vectors were calculated by computing their first

L = 100

entries. The error bars correspond to one standard deviation taken over ten repetitions of the algorithm. Runtimes correspond to a Python implementation of the algorithms running on 20 Intel Xeon CPUs operating at 2.4GHz with access to 252GB of shared memory.

These results show that for problems with the appropriate properties – namely very low rank and condition number – quantum-inspired algorithms can provide good estimates even for extremely high-dimensional problems. Importantly, these results are obtained in considerably less time than would be expected from the theoretical complexity bound $~ O (κ^{16} k^{6} ∥ A ∥_{F}^{6} / ε^{6})$ . Even setting $∥ A ∥_{F} = 1$ and $ε = 1$ , for $k = κ = κ_{β} = 10$ , the bound would suggest that an order of $10^{22}$ operations would be needed; this is not compatible with the runtimes we experience. This suggests that these theoretical bounds do not in fact reflect the practical runtime of the algorithms and therefore should not be used to evaluate their performance. Nevertheless, our findings do reflect that the quality of the estimates worsens considerably as rank and condition number are increased.

iii.2 Random matrices

We study the quantum-inspired algorithm for linear systems of equations $A x = b$ for randomly chosen $A$ and $b$

. To generate a random matrix

A

of dimension

m \times n

, rank

k

, and condition number

κ

, instead of creating

A

directly, we build the components of its singular value decomposition in matrix form,

A = U Σ V

. Here,

U

is an

m \times k

matrix whose columns are the left singular vectors of

A

. Similarly,

V

is a

k \times n

matrix whose rows are the right singular vectors, and

Σ

is a

k \times k

diagonal matrix whose entries are the corresponding singular values of

A

, which has rank

k

by construction. The details of numerical methods used to generate matrices

U

V

and

Σ

are explained in Appendix C. The vector

b

has been chosen to have the form

b = \sum_{ℓ = 1}^{k} β_{ℓ} u^{(ℓ)}

, with coefficients

β_{ℓ}

drawn from the standard normal distribution

N (0, 1)

We consider matrices consisting of $m = 40, 000$ rows and $n = 20, 000$ columns. These numbers are chosen to deal with the general case of rectangular matrices, while working at the limit of dimensions for which it would still be possible to calculate SVDs in a reasonable amount of time. This means, for simplicity, that in this setting we are not required to use specialized data structures to perform length-square sampling: we can instead store the entire distribution in memory and use fast built-in numerical sampling algorithms in Python.
In Figs. 2(a)-2(c) we monitor, respectively, the mean relative error $η_{σ} = \frac{1}{k} \sum_{ℓ = 1}^{k} \frac{| {~ σ}_{ℓ} - σ_{ℓ} |}{σ_{ℓ}}$ in estimating the singular values as well as the errors $η_{A} = ∥ ~ A - A ∥_{F} / ∥ A ∥_{F}$ and $η_{A}^{+} = ∥ {~ A}^{+} - A^{+} ∥_{F} / ∥ A^{+} ∥_{F}$ of the reconstructed matrices $A$ and $A^{+}$ . The approximation is accurate across all singular values even for a relatively small number of sampled rows and columns.

Figure 2: Errors in the approximate SVD of the random matrix $A$ of dimension $40, 000 \times 20, 000$ , rank $k = 5$ and condition number $κ = 5$ as the number of sampled rows $r$ and columns $c = r$ is increased. The figures show the error of the (a) singular values $η_{σ} = \sum_{ℓ = 1}^{k} | {~ σ}_{ℓ} - σ_{ℓ} | / σ_{ℓ}$ (b) reconstructed matrix $η_{A} = ∥ ~ A - A ∥_{F} / ∥ A ∥_{F}$ with $~ A = \sum_{ℓ = 1}^{k} {~ σ}_{ℓ} {~ u}^{(ℓ)} {~ v}^{{(ℓ)}^{T}}$ and (c) reconstructed pseudo-inverse $η_{A}^{+} = ∥ {~ A}^{+} - A^{+} ∥_{F} / ∥ A^{+} ∥_{F}$ with ${~ A}^{+} = \sum_{ℓ = 1}^{k} 1 / {~ σ}_{ℓ} {~ v}^{(ℓ)} {~ u}^{{(ℓ)}^{T}}$ . In all cases, error bars denote the standard deviation for 10 repetitions of the algorithm.

Empirically, the error in the approximation scales roughly as $1 / \sqrt{r}$ for these random matrices, which is in agreement with the matrix Chernhoff bound appearing in Theorem 3 of Ref. Gilyén et al. (2018). We also observe that the reconstructed matrix $~ A$ has smaller errors compared to the approximated pseudo-inverse matrix $A^{+}$ . This is because in reconstructing matrix $~ A$ , the largest singular values – which carry the dominant contributions – are approximated better by the FKV algorithm.

The main results are shown in Fig. 3, where we plot the error $η_{x}$ of the approximated solution vector $~ x$ calculated as the median of the relative errors $| {~ x}_{i} - x_{i} | / | x_{i} |$ with respect to all $n$ entries of vector $~ x$ . First, we investigate in Fig. 3(a) how the error $η_{x}$ changes as we increase the number of sampled rows $r$ and columns $c = r$ . All results have been obtained by taking $N = 10^{4}$ samples to estimate the coefficients $λ_{ℓ}$ defined in Eq. (12). As expected, $η_{x}$ decreases as the number of sampled rows and columns is increased. We have also superimposed error bars to show the standard deviation over 10 independent repetitions of the FKV algorithm.

Figure 3: Error $η_{x}$ of the solution vector for a $40, 000 \times 20, 000$ random matrix. (a) Error of the solution vector as a function of the number of sampled rows $r$ and columns $c$ entering the FKV algorithm, where we have fixed $c = r$ . We set a rank $k = 5$ and condition number $κ = 5$ . (b) Error in the solution vector as a function of matrix rank. There is a clear linear dependency as the matrix rank is increased. (c) Error in the solution vector as a function of the condition number of the input matrix. There is also a clear linear dependency, even for values spanning different orders of magnitude. In all cases, the coefficients $λ_{ℓ}$ have been estimated by performing $N = 10^{4}$ samples. Error bars show the standard deviation over 10 independent repetitions of the FKV algorithm.

In Fig. 3(b) we set $r = c = 4250$ and $κ = 5$ and investigate the error of the solution vector as a function of the rank $k$ of matrix $A$ . This plot reveals that the error $η_{x}$ increases linearly with rank. For $k = 5$ , a remarkably small relative error of $8.7 %$ is found for the entries of the approximated vector. However, we observe that the accuracy of the algorithm deteriorates rapidly for larger values of the rank. A similar trend is observed in Fig. 3(c) for the solution vector error as a function of the condition number $κ$ . Small errors of the order of $10 %$ are only obtained for small condition numbers, while large values of $η_{x} \geq 100 %$ already appear for $κ \geq 100$ . Note that this is the behaviour we would expect from Eq. (21) if the dominant source of errors originates from coefficient estimation. Similar plots for the error in the approximate SVD of matrix $A$ and in coefficients estimation are shown in Appendix C.

Finally, Tables 2 and 3 respectively summarize the errors and running time of the quantum-inspired algorithm for the case where the matrix $A$ has rank $k = 5$ and condition number $κ = 5$ . In Table 2, we show the average relative errors of singular values $η_{σ}$ and coefficients $η_{λ}$ as defined in the previous section. We also report the relative errors $η_{A} = ∥ A - ˜ A ∥_{F} / ∥ A ∥_{F}$ and $η_{A^{+}} = ∥ A^{+} - {˜ A}^{+} ∥_{F} / ∥ A^{+} ∥_{F}$ of the reconstructed matrices $˜ A$ and ${˜ A}^{+}$ . Sampling $4250$ rows and columns from the original $40000$ $\times$ $20000$ matrix $A$ is enough to have a good approximation full matrix. The largest error in Table 2 corresponds to the estimated coefficients.

Table 3 shows that the time to compute the approximate SVD using the FKV algorithm is orders-of-magnitude smaller than a direct calculation. Most of the overhead in running the quantum-inspired algorithms is associated with constructing the length-square distributions and estimating the coefficients. However, it is important to keep in mind that length-square distributions can in principle be constructed on-the-fly as the data from the matrix is generated, so it is best to separate this from the runtime of the actual algorithm. For this particular case of a matrix $A$ with low rank and small condition number, we find that the quantum-inspired algorithm outperforms the direct calculation.

	Parameters			Error
Case study	$r$	$c$	$N$	$η_{σ}$	$η_{A}$	$η_{A^{+}}$	$η_{λ}$	$η_{x}$
Random matrix	$4250$	$4250$	$10^{4}$	$0.010 \pm 0.005$	$0.028 \pm 0.004$	$0.101 \pm 0.027$	$0.387 \pm 0.191$	$0.087 \pm 0.053$

Table 2: Relative errors with corresponding standard deviations computed over

10

repetitions of the FKV algorithm of the quantum-inspired algorithm for linear systems

A x = b

A

is a random matrix with dimension

m = 40000

n = 20000

, rank

k = 5

, and condition number

κ = 5

. Approximated singular vales and vectors were obtained by decomposing matrix

C

with dimension

r = 4250

and

c = 4250

. The coefficients

˜ λ

have been estimated by performing

N = 10^{4}

samples. Here,

η_{σ}

is the error in approximating singular values

~ σ

, while

η_{A}

and

η_{A^{+}}

are the errors of the reconstructed matrices

˜ A

and

{˜ A}^{+}

, respectively.

η_{λ}

is the error in estimating coefficients, and

η_{x}

is the the median of the relative errors

| {~ x}_{i} - x_{i} | / | x_{i} |

with respect to all

n = 20000

entries of the approximate vector

~ x

	Quantum-inspired algorithm					Direct calculation
Case study	$t_{L S}$	$t_{S V D}^{C}$	$t_{λ}$	$t_{x}$	$t_{t o t a l}$	$t_{S V D}^{A}$	$t_{λ}$	$t_{x}$	$t_{t o t a l}$
Random matrix	$1488.8$	$83.9$	$554.7$	$343$	$2470.4$	$5191.1$	$1.4$	$0.0003$	$5192.5$

Table 3: Running times (in seconds) for the quantum-inspired algorithm applied to solve the system of linear equations

A x = b

with the same parameters as reported in Table 2. The parameter

t_{L S}

is the time required to construct length-squared (LS) probability distributions over rows and columns.

t_{S V D}^{C}

accounts for the time spent sampling rows and columns from matrix

A

and building and decomposing the sampled matrix

C

t_{λ}

is the running time of the estimation of the coefficients

λ

. Analogous quantities are defined for the direct calculation method.

t_{x}

is the elapsed time to sample

500

entries of the approximated solution vector in the quantum-inspired algorithm and the time spent to compute the exact solution vector in the direct calculation method. Runtimes correspond to a Python implementation of the algorithms running on two Intel Xeon CPUs operating at 2.4GHz with access to 252GB of shared memory.

iii.3 Portfolio optimization

We study an application of the quantum-inspired algorithm for linear systems of equations to a canonical financial problem: portfolio optimization. In its simplest version, commonly known as the Markowitz mean-variance model Markowitz (1952), the goal is to optimally invest wealth across $n$ possible assets that are modeled only by their expected returns and correlations. Consider a vector of returns $r_{j} = (r_{1, j}, r_{2, j}, \dots, r_{n, j})^{T}$ , where $r_{i, j}$ is the return of asset $i$ on the $j$ -th day. The vector of expected returns $r$ and the correlation matrix $Σ$ are defined as

r = \frac{1}{n} n \sum j = 1 r_{j}, Σ = \frac{1}{n} n \sum j = 1 r_{j} r_{j}^{T} .

(30)

The vector $r$ captures the expected returns for each asset, while the correlation matrix $Σ$ expresses fluctuations around mean values, which in the model is interpreted as risk. An optimal portfolio in this setting is one that minimizes risk for any given target expected return. This corresponds to the optimization problem

	$Minimize: w w^{T} Σ w$		(31)
	$subject to: r^{T} w = μ,$		(32)

where $w$ is a portfolio allocation vector that determines how much wealth is invested in each asset. We allow $w$ to have negative entries, which is interpreted as an indication to short-sell the corresponding asset. As discussed in Ref. Rebentrost and Lloyd (2018), this problem can be equivalently cast as a linear system of equations:

[\begin{matrix} 0 & r^{T} r & Σ \end{matrix}] [\begin{matrix} ν w \end{matrix}] = [\begin{matrix} μ 0 \end{matrix}],

(33)

where $ν$ is a Lagrange multiplier. In this case, the linear system $A x = b$ is given by $A = [\begin{matrix} 0 & r^{T} r & Σ \end{matrix}]$ , $b = [\begin{matrix} μ 0 \end{matrix}]$ , and $x = [\begin{matrix} ν w \end{matrix}]$ . To compute a solution to the linear system, we employ the Moore-Penrose pseudo-inverse $A^{+}$ , which can be expressed in terms of the singular values and vectors of $A$ as $A^{+} = \sum_{ℓ = 1}^{k} \frac{1}{σ_{ℓ}} v^{(ℓ)} {u^{(ℓ)}}^{T}$ . The solution vector $x$ is then

	$x$	$= A^{+} b = k \sum ℓ = 1 \frac{1}{σ_{ℓ}} ⟨ u^{(ℓ)}, b ⟩ v^{(ℓ)}$
		$= k \sum ℓ = 1 \frac{1}{σ_{ℓ}^{2}} ⟨ v^{(ℓ)}, A^{T} b ⟩ v^{(ℓ)} = k \sum ℓ = 1 λ_{ℓ} v^{(ℓ)},$		(34)

where we have identified the coefficients $λ_{ℓ} = ⟨ v^{(ℓ)}, A^{T} b ⟩ / σ_{ℓ}^{2}$ as in Eq. (3). The vector $b$ has only one non-zero entry, so it holds that $A^{T} b = μ A_{1}$ , where $A_{1}$ is the first row of $A$ .

To test the algorithm in this setting, we employ publicly-available pricing data for the stocks comprising the S&P 500 stock index during the five-year period 2013-2018 Nugent (2018). Since not all stocks remain in the index during this time, we restrict to the 474 stocks that do, leading to a matrix $A$ of dimension $475 \times 475$ . Returns are calculated on a daily basis based on opening price. The matrix has full rank of $k = 475$ , its largest and smallest singular values are respectively $σ_{m a x} = 40.2$ , $σ_{m i n} = 1.8 \times 10^{- 4}$ , and its condition number is $κ = 2.23 \times 10^{4}$ . More details can be found in Appendix D. We set a target return $μ$ equal to the average return over all 474 stocks in the index.

In Tables 4 we report the errors of the quantum-inspired algorithm. The SVD of matrix $A$ was approximated by sampling $340$ rows and columns from the full matrix. As in the previous example, we took $N = 10^{4}$ samples to estimate the coefficients $λ_{ℓ}$ . The relative errors characterizing the approximate solution vector are considerable, with values of $η_{x} = 0.74$ . In this example, we have calculated the approximate solution vector by using Eq. (III.3) within a low-rank approximation with $k = 10$ . Notice from Table 4 that the large errors $η_{A}^{+}$ in the reconstructed pseudo-inverse matrix and estimated coefficients $η_{λ}$ translate into large relative errors of the solution vector.

Finally, in Table 5 we report the running-times for the algorithm. Since the matrices involved in this example are small, exact and approximate SVDs can be computed quickly, but for the quantum-inspired algorithm, there is significant overhead due to coefficient estimation and sampling from the solution vector.

	Parameters			Error
Case study	$r$	$c$	$N$	$η_{σ}$	$η_{A}$	$η_{A}^{+}$	$η_{λ}$	$η_{x}$
S&P 500 stock index	$340$	$340$	$10^{4}$	$0.08 \pm 0.02$	$0.16 \pm 0.03$	$1.13 \pm 0.12$	$1.58 \pm$

Quantum-inspired algorithms in practice

Authors

Quantum-Inspired Classical Algorithms for Singular Value Transformation

Parametrized Complexity of Quantum Inspired Algorithms

Quantum-inspired sublinear classical algorithms for solving low-rank linear systems

A quantum-inspired classical algorithm for recommendation systems

On exploring practical potentials of quantum auto-encoder with advantages

Quantum-Inspired Algorithms from Randomized Numerical Linear Algebra

Classical Algorithms from Quantum and Arthur-Merlin Communication Protocols

Code Repositories