Robust Dequantization of the
Quantum Singular Value Transformation and
Quantum Machine Learning Algorithms

Our first contribution is to develop a framework for working with vectors and matrices given by approximate sampling-and-query access. Following the approach developed in [5] (which used the concept of oversampling-and-query access), we develop our framework by introducing the concept of approximate over-sampling-and-query access. In this setting, instead of being able to generate samples from a distribution $\tilde{p}_{u}$ close to $p_{u}$, we can only generate samples from a distribution $\tilde{q}_{u}\colon\{1,\ldots,n\}\to[0,1]$ such that $\tilde{q}_{u}(i)\geq\frac{1}{\phi}\tilde{p}_{u}(i)$ holds for all $i\in\{1,\ldots,n\}$, for some $\phi\geq 1$ called the oversampling parameter. (Note that taking $\phi=1$ gives the usual approximate sampling-and-query access.) We prove the following two properties for approximate oversampling-and-query access.

• •

  Closeness property: given approximate oversampling-and-query access to a small number of vectors, we can implement efficiently approximate oversampling-and-query access to a linear combination of these vectors (with a slightly worse oversampling parameter).

• •

  Conversion property: given approximate oversampling-and-query access to a vector, we can implement efficiently approximate sampling-and-query access with an overhead in the complexity corresponding to the oversampling parameter.

These two properties generalize the two properties shown (for the setting of perfect sampling-and-query access) in [5]. The closeness property, which usually does not hold for approximate sampling-and-query access, is the reason why the notion of approximate oversampling-and-query access is very convenient to deal with. In addition, the conversion property shows that working with oversampling-and-query access is meaningful, since at the end of the computation it is possible to recover sampling-and-query access from the output (with an overhead in the complexity corresponding to the oversampling parameter).

Our second contribution is to develop a robust version of one of the central techniques in the area of randomized numerical linear algebra: matrix multiplication based on importance matrix sketches. Given two matrices $X\in\mathbb{C}^{m\times n}$ and $Y\in\mathbb{C}^{m\times n^{\prime}}$, the approach by Drineas, Kannan and Mahoney [13], on which the dequantization framework by Chia, Gilyén, Li, Lin, Tang and Wang [5] is also based, creates a matrix $\Sigma\in\mathbb{R}^{r\times m}$ with only one nonzero entry per row, and estimates the matrix product $X^{\dagger}Y$ by the matrix product $X^{\dagger}\Sigma^{\dagger}\Sigma Y$. Note that $X^{\dagger}\Sigma^{\dagger}$ is an $n\times r$ matrix obtained by selecting $r$ columns of $X^{\dagger}$ (i.e., $r$ rows of $X$) and renormalizing them, while $\Sigma Y$ is an $r\times n^{\prime}$ matrix obtained by selecting $r$ rows of $Y$ and renormalizing them. The matrix product $X^{\dagger}\Sigma^{\dagger}\Sigma Y$ can thus be considered as a “sketch” of the matrix product $X^{\dagger}Y$. The quality of the approximation, i.e., the Frobenius norm

$$\left\|X^{\dagger}\Sigma^{\dagger}\Sigma Y-X^{\dagger}Y\right\|_{\mathrm{F}}$$

depends naturally on the choice of $\Sigma$. Drineas, Kannan and Mahoney [13] analyzed the optimal choice for $\Sigma$, and showed that this $\Sigma$ can be constructed using (length-squared) sampling-and-query access to the matrices $X$ and $Y$.⁵⁵5Sampling-and-query access to a matrix is defined using the concept of sampling-and-query access to a vector: we say that we have sampling-and-query access to a matrix $A\in\mathbb{C}^{m\times n}$ if we have sampling-and-query access to each row of $A$ and we have sampling-and-query access to the vector of row norms of the matrices $A$, i.e., the vector $(\left\|A(1,\cdot)\right\|,\ldots,\left\|A(m,\cdot)\right\|)$. This notion has been extended to oversampling-and-query access to a matrix in the literature [5, 13]. In this work we further extend it to $\varepsilon$-approximate (over)sampling-and-query access to a matrix (see Definitions 7 and 8 in Section 2.3). In the dequantization setting, this result has been used in particular to obtain a sketch of the singular value transformation of a matrix [5, Theorem 5.1].

In our setting, where we only have approximate (over)sampling-and-query access to the matrices $X$ and $Y$, several challenges appear when trying to adapt this technique. First, constructing efficiently unbiased estimators for matrix multiplication seems hopeless since the probability distribution from which we can sample is biased. Second, and more importantly, it is even unclear whether the approach from [5, 13] can be adapted, for the same reasons as mentioned in Section 1.2 for the inner product (observe that inner product is a special case of matrix product).

Our main technical contribution is to show that a robust version of matrix multiplication is possible, with a small bias depending on the parameter $\varepsilon$. Here are informal versions of our two results related to matrix multiplication (we refer to Sections 4.1 and 4.2 for the formal statements):

Theorem 2 gives a smaller bias and an exponentially better dependence on $\delta$ than Theorem 1. On the other hand, it requires $\varepsilon$-approximate sampling-and-query access to both $X$ and $Y$. We note that the dependence in $\delta$ in Theorem 1 is likely to be inevitable since even with perfect sampling-and-quantum access, it is not known how to achieve a better dependence on $\delta$ when given access to only one matrix [13].

We then use these techniques to show how to compute a sketch of the singular value transformation of a matrix, similarly to what [5] has done in the perfect sampling-and-query access setting. Here is the informal version of our result (see Section 4.3 for the formal version):

Theorem 3 shows that in order to get a good “sketch” of a singular value transformation $f(X^{\dagger}X)$ of the matrix $X^{\dagger}X$ it is enough to compute $\bar{f}(CC^{\dagger})$ for the much smaller matrix $CC^{\dagger}\in\mathbb{C}^{r\times r}$. This technique is the basis for most of our advanced applications.

By applying all the above techniques, we are able to robustly dequantize (i.e., dequantize even when the input is given via $\varepsilon$-approximate sampling-and-query access for $\varepsilon>0$) most⁶⁶6Due to space constraints, in this paper we only present our robust dequantization results for a few important quantum algorithms. Our framework does apply to more quantum algorithms. Actually, we have not been able to find an example of quantum algorithm that can be dequantized but not robustly dequantized. of the quantum algorithms that have been dequantized in the literature. Below we briefly discuss all our applications and refer to Section 6 for details. Simplified statements of the complexities of our results and their comparison with prior works can be found in Table 1.⁷⁷7In this paper, the notation $\tilde{O}(\cdot)$ removes factors polylogarithmic in the dimension of the matrices and vectors.

• •

  Inner product (Section 3 and Section 6.1). In the quantum setting, the SWAP test can be used to estimate the absolute value of the inner product between two pure states with additive error $\xi$ using $O(1/\xi^{2})$ copies of the states.⁸⁸8If quantum circuits constructing the two quantum states are available, then the complexity can be further reduced to $O(1/\xi)$ using amplitude estimation. The dequantization technique from [5, 39] estimates the inner product of two vectors $u$ and $v$ with error $\xi\left\|u\right\|\left\|v\right\|$ in time $O(1/\xi^{2})$ given sampling-and-query access to one of the vectors and query access to the other, which matches the complexity of the SWAP test for unit-norm vectors. By applying Theorem 1, we immediately obtain a classical algorithm that estimates the inner product with error $(2\sqrt{2\varepsilon}+\xi)\left\|u\right\|\left\|v\right\|$ in time $O(1/\xi^{2})$ when given $\varepsilon$-approximate sampling-and-query access to one of the vectors and query access to the other (Proposition 1 in Section 3 and Proposition 6 in Section 6.1). The term $2\sqrt{2\varepsilon}$ here is a bias (coming from the bias in Theorem 1) due to the fact that we are sampling from a biased distribution.⁹⁹9Such a bias is inevitable when estimating the inner product using samples from a biased distribution. Note that a similar bias would appear for quantum algorithms based on the SWAP test in the analogue setting where the input states are at (trace) distance $\varepsilon$ from the ideal states. This result shows that when $\varepsilon$ is small enough, i.e., when we are sampling from a distribution close enough to the ideal distribution in the total variation distance, we can obtain a performance close to the performance of algorithms having access to the ideal distribution.

• •

  “Low-rank” Quantum Singular Value transformation (Section 6.2). Chia, Gilyén, Li, Lin, Tang and Wang [5] have shown how to dequantize in the low-rank regime the Quantum Singular Value Transformation (QSVT), which is a recent powerful paradigm developed by Gilyén, Su, Low and Wiebe [21] that can be used to recover most known quantum algorithms and construct new ones (we refer to [33] for a good survey). Given an even polynomial $\mathfrak{p}$ of degree $d$, a matrix $A\in\mathbb{C}^{m\times n}$ and a vector $b$ with $\left\|b\right\|=1$, the QSVT is a technique to construct a good approximation of a quantum state proportional to the vector $\mathfrak{p}(\sqrt{A^{\dagger}A})b$. In the low-rank regime, which corresponds to the assumption $\left\|A\right\|_{\mathrm{F}}=1$ (or more generally $\left\|A\right\|_{\mathrm{F}}\leq 1$), the complexity of the approach is $\tilde{O}(d/\|\mathfrak{p}(\sqrt{A^{\dagger}A})b\|)$. Chia, Gilyén, Li, Lin, Tang and Wang [5] showed how to dequantize this technique by proving that given sampling-and-query access to $A$ and $b$, it is possible to implement with high probability sampling-and-query access to a vector close to $\mathfrak{p}(\sqrt{A^{\dagger}A})b$ in time $\tilde{O}(d^{16}/\|\mathfrak{p}(\sqrt{A^{\dagger}A})b\|^{6})$. By adapting the methodology from [5], we show how to achieve the same running time when given only $\varepsilon$-approximate sampling-and-query access to $A$ and $b$ for $\varepsilon>0$, when $\varepsilon$ is small enough.

• •

  Sparse Quantum Singular Value Transformation (Section 6.3). Gharibian and Le Gall [18] have considered the problem of estimating the value $v^{\dagger}\mathfrak{p}(\sqrt{A^{\dagger}A})u$ for a matrix $A\in\mathbb{C}^{m\times n}$ such that $\left\|A\right\|=1$ (or more generally $\left\|A\right\|\leq 1$) and two vectors $u,v\in\mathbb{C}^{n}$. Note that the condition $\left\|A\right\|\leq 1$ is significantly weaker than the condition $\left\|A\right\|_{\mathrm{F}}\leq 1$ since the inequality $\left\|A\right\|\leq\left\|A\right\|_{\mathrm{F}}$ always holds but $\left\|A\right\|_{\mathrm{F}}$ can be as large as $\sqrt{n}\left\|A\right\|$. On the other hand, with this weaker assumption the problem can only be solved efficiently in the classical setting for special types of matrices, such as $s$-sparse matrices (i.e., matrices containing at most $s$ nonzero entries per row and column). The main result from [18] indeed shows that when given sampling-and-query access to an $s$-sparse matrix $A$, a good estimate of the value $v^{\dagger}\mathfrak{p}(\sqrt{A^{\dagger}A})u$ can be computed in $\tilde{O}(s^{d})$ time, which matches (up to polynomial factors) the complexity of quantum algorithms based on the QSVT for constant $s$ and constant $d$ (Ref. [18] also shows that the problem becomes $\mathsf{BQP}$-hard for $d=\mathrm{poly}(n)$ and $s$ constant). By combining our techniques for robust estimation of the inner product with the main technical result from [18], we show that this dequantization result even holds when considering $\mathsf{SQ}_{\varepsilon}(v)$ for $\varepsilon>0$, when $\varepsilon$ is small enough.

  Using this result, we can generalize the other dequantization results in [18, Section 4], and in particular obtain an efficient classical algorithm computing a constant-precision estimation of the ground energy of a local Hamiltonians when given $\varepsilon$-approximate sampling-and-query access to a guiding vector $u$ that has constant overlap with the ground state (Theorem 1 in [18] showed this result only for the case $\varepsilon=0$). Similarly, the implications to the Quantum PCP conjecture and the No Low-Energy Samplable States (NLSS) conjecture discussed in [18, Section 5] also generalize to the case of approximate sampling-and-query access to the witness states.

  While not considered in [18], the same approach can be used for sparse matrix inversion as well. As discussed in [21] (see also [33]), the QSVT can be used for matrix inversion by choosing a polynomial $\mathfrak{p}$ that approximates the function $f(x)=1/x$ on the interval $[-1,1]\setminus(-1/\kappa,1/\kappa)$, where $\kappa$ is the condition number. The analysis from [21] shows that a good enough approximation can be obtained by a polynomial of degree $d=O(\kappa\log(\kappa/\eta))$, where $\eta$ is the needed precision. The dequantization technique from [18] then leads to a classical algorithm with complexity $\tilde{O}(s^{\kappa\log(\kappa/\eta)})$ given exact sampling-and-query access to the input. Our result shows that this complexity holds even when considering approximate sampling-and-query access.

• •

  Supervised clustering (Section 6.4). Lloyd, Mohseni, and Rebentrost [31] developed a quantum algorithm for the supervised clustering problem, an important problem in machine learning. From a computational perspective, this quantum algorithm simply computes the inner product between two vectors representing appropriately the data, which can be done in $\tilde{O}\left(1/\eta\right)$ time using the SWAP test, where $\eta$ represents the precision parameter. Prior dequantization works [5, 40] showed how to construct classical algorithms solving this problem in time $O\left(1/\eta^{2}\right)$ given sampling-and-query access to the input. By applying our algorithm for robust estimation of the inner product, we immediately obtain the same complexity when given only $\varepsilon$-approximate sampling-and-query access for $\varepsilon>0$.

• •

  Recommendation systems (Section 6.5). In recommendation systems, the main computational task reduces to sampling from a row of a low-rank matrix close to the input matrix $A$ (which corresponds to the users’ preference matrix). A key parameter for expressing the complexity of quantum algorithms for recommendation systems [29] and their dequantization [5, 39] is the ratio $K=\left\|A\right\|_{\mathrm{F}}^{2}/\sigma^{2}$, where $\sigma$ is a threshold value used for the specification of the problem. Kerenidis and Prakash [29] have shown how to solve the problem in time $\tilde{O}(\sqrt{K})$ on a quantum computer when $A$ is given in QRAM. Ref. [5] gave a classical algorithm solving the problem in time $\tilde{O}(K^{16})$ when $A$ is given via sampling-and-query access, which matches the complexity of the quantum algorithm up to a (large) polynomial factor and improved the complexity of Tang’s original algorithm [39]. By applying our framework, we show how to obtain the same complexity when given only $\varepsilon$-approximate sampling-and-query access to the input for $\varepsilon>0$, when $\varepsilon$ is small enough.

• •

  “Low-rank” matrix inversion (Section 6.6). While the matrix inversion problem solved by the HHL quantum algorithm is $\mathsf{BQP}$-hard [23] and thus unlikely to be dequantized, several low-rank versions of matrix inversion (or, more precisely, solving a system of linear equations) have been studied in recent works on dequantization [4, 5, 22]. Ref. [5], in particular, showed that classical algorithms can solve linear equations in time $\tilde{O}(K^{14})$ given sampling-and-query access to the matrix and the vector specifying the system of linear equations. Since quantum algorithms solve this problem in time $\tilde{O}(\sqrt{K})$ when $A$ is given in QRAM [21, 22], the classical algorithms match this complexity up to a (large) polynomial factor. By applying our framework, we show how to obtain the same complexity when given only $\varepsilon$-approximate sampling-and-query access for $\varepsilon>0$, when $\varepsilon$ is small enough.

All these results give strong evidence for the lack of exponential speedup for quantum machine learning algorithms in the QRAM model (and additionally for the framework from [18]) even when the classical algorithms only have approximate sampling-and-query access to the input. The main conceptual message of this paper is thus that no quantum advantage seems to arise in machine learning (at least for the problems we consider) from the ability of quantum algorithms to natively robustly manipulate classical information encoded as quantum states.

Very recently, Bakshi and Tang [2] have constructed quantum-inspired classical algorithms for the Quantum Singular Value transformation and its applications with improved bounds compared with [5]. For low-rank matrix inversion, recent works by Gilyén, Song and Tang [20] and Montanaro and Shao [38] have given quantum-inspired classical algorithms with improved bounds over [5] as well. These improved algorithms, which use iterative methods and new classical techniques, assume perfect sampling-and-query access to the data, as in [5]. It would be interesting to obtain similar improvements in the case of approximate sampling-and-query access.

We describe the main idea behind our robust estimation of the inner product. Consider two unit-norm vectors $u,v\in\mathbb{R}^{n}$, where $u$ is given via $\varepsilon$-approximate sampling-and-query access and $v$ is given by query access (i.e., on an input $i\in\{1,\ldots,n\}$ we can query $v(i)$). A straightforward adaptation of the method described in Section 1.2 would be: sample an index $i\in\{1,\ldots,n\}$ according to the distribution $\tilde{p}_{u}$ and output the value $v(i)/u(i)$. Comparing with Equation (1), we see that indices $i$ such that $u(i)$ is very small but $\tilde{p}_{u}(i)$ is large lead to large bias. This suggests the followings strategy: use instead an estimator of the type

$$\begin{cases}0&\textrm{ if }i\in\Gamma,\\ \frac{v(i)}{u(i)}&\textrm{ if }i\notin\Gamma,\end{cases}$$

for a set $\Gamma=\{i\in\{1,\ldots,n\}\>|\>\left|u(i)\right|\leq\theta(i)\}$, where $\theta(i)$ is a well chosen threshold function. Ideally, we would like to choose $\theta(i)$ based on the value of $\tilde{p}_{u}(i)$. This is unfortunately not possible since we do not have access to the value $\tilde{p}_{u}(i)$. Our key discovery is that choosing $\theta(i)$ based on the value of $v(i)$ does work. Taking such as estimate leads to a biased estimator where the bias comes from two parts: the first part coming from the fact that we are sampling from $\tilde{p}_{u}$ instead of $p_{u}$ and the second part coming from the fact that we disregarding the contribution of the indices in the set $\Gamma$. When setting $\theta(i)=\gamma\left|v(i)\right|$ for some constant $\gamma$, we discover that it is possible to derive “complementary” upper bounds on these two biases: the upper bound on the first bias is proportional to $\gamma$ while the upper bound on the second bias is inverse-proportional to $\gamma$. We then take the optimal value for $\gamma$ and set $\theta(i)$ accordingly. Concretely, we take $\theta(i)=\sqrt{2\varepsilon}\left|v(i)\right|$ (with additional normalization factors when $u$ and $v$ are not unit-norm), which gives a bias of $\sqrt{2\varepsilon}+\sqrt{2\varepsilon}=2\sqrt{2\varepsilon}$. We refer to Section 3 for a fairly self-contained presentation of the details of our technique for robust estimation of the inner product.¹⁰¹⁰10The inner product estimation algorithm of Section 3 is actually a special case of Theorem 1 in Section 4.1. We present it in Section 3 as a “warm-up” since for this special case the presentation is significantly lighter (in particular, we do not need the notion of matrix sketch introduced in in Section 4.1).

Our second main technical contribution is extending this strategy to compute a general matrix product $XY$. We adapt the notion of matrix sketch, which is one of the central concepts in linear algebra algorithms based on random sampling [11, 12, 13, 14, 15, 17, 27, 28], to the setting of approximate sampling (see Definition 9 in Section 4.1). This is the key step that enables us to adapt the matrix multiplication sketching technique by Drineas, Kannan and Mahoney [13] to the setting where only $\varepsilon$-approximate sampling-and-query access to one matrix is given, for $\varepsilon>0$. When given approximate sampling access only to one of the matrices ($X$ or $Y$), we apply a technique similar to our robust evaluation of the inner product (indeed, we obtain the $2\sqrt{2\varepsilon}$ term in the bias of the sketch in Theorem 1). When given approximate sampling access to both matrices $X$ and $Y$, we adapt the technique based on joint-matrix sampling from [13], which was also used in [5]. The first step is again to adapt the definition to the setting of approximate sampling (Definition 12 in Section 4.2). Using the inequality of arithmetic and geometric means and a finer analysis of the success probability based on MacDiarmid’s bounded difference inequality, this technique leads to an algorithm (Theorem 2) with exponentially better dependence in the error parameter $\delta$.

Our third main technical contribution, crucial to dequantize quantum algorithms for the low-rank quantum singular value transformation and recommendation systems, is the conversion technique mentioned in Section 1.3, which converts approximate oversampling-and-query access into approximate sampling-and-query. A similar conversion has been used in the perfect setting in [5, 39], based on the standard technique in probabilistic algorithms called rejection sampling (see, e.g., [9] for a detailed description of this technique). Our main technical contribution is to develop a robust version of rejection sampling.

Standard rejection sampling is a method to generate samples from a (hard) distribution $p_{1}\colon\{1,\ldots,n\}\to[0,1]$ given the ability to generate samples from another (easier) distribution $p_{2}\colon\{1,\ldots,n\}\to[0,1]$. The method works as follows: sample $j\in\{1,\ldots,n\}$ from $p_{2}$; output it with probability $\frac{p_{1}(j)}{mp_{2}(j)}$ for a large value $m$ (large enough such that $p_{1}(j)\leq mp_{2}(j)$ holds for all $j\in\{1,\ldots,m\}$) and otherwise report “failure”. It is easy to see that the probability of outputting $j$ is $p_{1}(j)/m$ and the probability that the process does not fail is $1/m$. In consequence, conditioned on the event that the process does not fail, the probability that sample $j$ is output is precisely $p_{1}(j)$. Repeating the process $\Theta(m)$ times will then output a sample drawn from the distribution $p_{1}$ with high probability. In Section 5.1, we show how to generalize this technique to the case where we can only generate samples from a probability distribution $\tilde{p}_{2}$ that is close to $p_{2}$ in total variation distance. At first glance, it seems that the same difficulty as when trying to extend Equation (1) to the approximate sampling setting arises. We show that the condition $p_{1}(j)\leq mp_{2}(j)$ is actually strong enough to control the bias, which enables us to sample with high probability from a probability distribution close to $p_{1}$ in total variation distance (see Proposition 3 in Section 5.1).

In this paper we assume that arithmetic operations in $\mathbb{C}$ (i.e., addition, subtraction, multiplication and division of two scalars) can be done at unit cost.¹¹¹¹11All our upper bounds can easily be converted into upper bounds for the bit-complexity of the problem as well if the notions of query access and sampling access to the input are redefined appropriately to take in consideration the bit-complexity of the input vectors and matrices. For any complex number $x$, we write $x^{\ast}$ its complex conjugate. We denote by $\mathbb{R}_{\geq 0}$ the set of non-negative real numbers. We say that a polynomial $\mathfrak{p}\in\mathbb{C}[x]$ is even if $\mathfrak{p}(x)=\mathfrak{p}(-x)$ for all $x\in\mathbb{C}$. Given a function $f\colon\mathbb{R}\to\mathbb{C}$, some set $\Lambda\subseteq\mathbb{R}$ and a positive real number $L$, we say that $f$ is $L$-Lipschitz on $\Lambda$ if the inequality $\left|f(x)-f(y)\right|\leq L\left|x-y\right|$ holds for any $x,y\in\Lambda$.

For a matrix $A\in\mathbb{C}^{m\times n}$, we write its entries as $A(i,j)$ for $(i,j)\in\{1,\ldots,m\}\times\{1,\ldots,n\}$. For each $i\in\{1,\ldots,m\}$, we denote the $i$-th row of $A$ as $A(i,\cdot)$ . For any $j\in\{1,\ldots,n\}$, we denote the $j$-th column of $A$ as $A(\cdot,j)$. We write $A^{\ast}$ for its conjugate and $A^{\dagger}$ for its transpose-conjugate, i.e., the matrices such that $A^{\ast}(i,j)=A(i,j)^{*}$ and $A^{\dagger}(i,j)=A(j,i)^{*}$ for all $(i,j)\in\{1,\ldots,m\}\times\{1,\ldots,n\}$. We use similar notations for vectors.

For a vector $u\in\mathbb{C}^{n}$, we denote its $\ell_{2}$ norm as $\left\|u\right\|$. For a matrix $A\in\mathbb{C}^{m\times n}$, we denote its spectral norm as $\left\|A\right\|$, and denote its Frobenius norm as $\left\|A\right\|_{\mathrm{F}}$, i.e.,

$$\left\|A\right\|_{\mathrm{F}}=\sqrt{\sum_{i=1}^{m}\sum_{j=1}^{n}\left|A(i,j)\right|^{2}}.$$

We have $\left\|A\right\|\leq\left\|A\right\|_{\mathrm{F}}\leq\sqrt{\min(m,n)}\left\|A\right\|$.¹²¹²12In this paper we systematically use the upper bound $\left\|A\right\|\leq\left\|A\right\|_{\mathrm{F}}$ to state the bounds on the complexity of our algorithms only in term of $\left\|A\right\|_{\mathrm{F}}$. It is possible (and done in prior works [5]) to give a finer analysis of the statements of upper bounds using both $\left\|A\right\|$ and $\left\|A\right\|_{\mathrm{F}}$ but for clarity we abstain from doing this since this gives more complicated bounds.

Given a Hermitian matrix $H\in\mathbb{C}^{n\times n}$, we order its eigenvalues in non-increasing order and denote them $\lambda_{1}(H)\geq\lambda_{2}(H)\geq\cdots\geq\lambda_{n}(H)$. We write

$$\mathrm{Spec}(H)=\{\lambda_{1}(H),\lambda_{2}(H),\ldots,\lambda_{n}(H)\}$$

for the spectrum of $H$.

In Section 4.3 we will need the following two lemmas.

For a positive semidefinite Hermitian matrix $H\in\mathbb{C}^{n\times n}$, consider its spectral decomposition

$$H=\sum_{i=1}^{n}\lambda_{i}(H)v_{i}v_{i}^{\dagger}.$$

where the $\{v_{i}\}_{i}$ are orthonormal vectors in $\mathbb{C}^{n}$. For any function $f\colon\mathbb{R}_{\geq 0}\to\mathbb{C}$, the matrix

$$f(H)=\sum_{i=1}^{n}f(\lambda_{i}(H))v_{i}v_{i}^{\dagger}$$

is called the singular value transformation of $H$ associated to $f$ (which in this case coincides with the eigenvalue transformation of $H$ associated to $f$). We use this definition to define the concept, which we will consider in Sections 6.2 and 6.3, of (even) polynomial transformation of an arbitrary (i.e., non-Hermitian) matrix. For an even polynomial $\mathfrak{p}\in\mathbb{C}[x]$ and a matrix $A\in\mathbb{C}^{m\times n}$, the (even) polynomial transformation of $A$ associated to $\mathfrak{p}$ is defined as $\mathfrak{p}(\sqrt{A^{\dagger}A})$. Let us write

$$\mathfrak{p}(x)=a_{0}+a_{2}x^{2}+a_{4}x^{4}+\cdots+a_{d}x^{d},$$

where $d$ is the (even) degree of $\mathfrak{p}$ and $a_{0},\ldots,a_{d}$ are coefficients. It is easy to check that the following relation holds:

$$\mathfrak{p}(\sqrt{A^{\dagger}A})=a_{0}I_{n}+a_{2}A^{\dagger}A+a_{4}(A^{\dagger}A)^{2}+\cdots+a_{d}(A^{\dagger}A)^{d},$$

where $I_{n}$ denotes the identity matrix of dimension $n$.

We now define the singular value transformation of an arbitrary matrix. This concept will be used in Sections 6.5 and 6.6. For any matrix $A\in\mathbb{C}^{m\times n}$, consider its spectral singular value transformation

$$A=\sum_{i=1}^{\min(m,n)}\sigma_{i}u_{i}v_{i}^{\dagger},$$

where the $\sigma_{i}$’s are nonnegative real numbers, $\{u_{i}\}_{i}$ are orthonormal vectors in $\mathbb{C}^{m}$ and $\{v_{i}\}_{i}$ are orthonormal vectors in $\mathbb{C}^{n}$. For any function $f\colon\mathbb{R}_{\geq 0}\to\mathbb{C}$ such that $f(0)=0$, the matrix

$$f(A)=\sum_{i=1}^{\min(m,n)}f(\sigma_{i})u_{i}v_{i}^{\dagger}$$

is called the singular value transformation of $A$ associated to $f$. Note that when $A$ is a positive semidefinite Hermitian matrix we recover the above definition.¹³¹³13For arbitrary matrices, the condition $f(0)=0$ guarantees that the singular value transformation is well defined.

Given a probability distribution $p\colon\{1,\ldots,n\}\to[0,1]$, we denote its support as $\mathrm{supp}(p)$, i.e.,

$$\mathrm{supp}(p)=\{i\in\{1,\ldots,n\}\>|\>p(i)>0\}.$$

Given two probability distributions $p,q\colon\{1,\ldots,n\}\to[0,1]$ we write their total variation distance as

$$|p-q|_{\mathrm{tv}}=\frac{1}{2}\sum_{i=1}^{n}\left|p(i)-q(i)\right|.$$

As already mentioned in the introduction, for any nonzero vector $u\in\mathbb{C}^{n}$ we define the probability distribution $p_{u}\colon\{1,\ldots,n\}\to[0,1]$ as

$$p_{u}(i)=\frac{|u(i)|^{2}}{\left\|u\right\|^{2}}$$

for each $i\in\{1,\ldots,n\}$.

We will often work with matrix-valued random variables. We recall the definition of the expectation and variance of such variables. A matrix-valued random variable $Z\in\mathbb{C}^{m\times n}$ is a matrix where each entry $Z(i,j)$ is a complex-valued random variable, i.e., $Z(i,j)=X(i,j)+i\>Y(i,j)$ for two real-valued random variables $X(i,j)$ and $Y(i,j)$. The expectation of $Z$ is the matrix $C\in\mathbb{C}^{m\times n}$ such that the $C(i,j)=\mathbb{E}\!\left[X(i,j)\right]+i\>\mathbb{E}\!\left[Y(i,j)\right]$. The variance of $Z$ is

$$\mathrm{Var}\left[Z\right]=\mathbb{E}\!\left[\left\|Z-\mathbb{E}\!\left[Z\right]\right\|_{\mathrm{F}}^{2}\right]=\mathbb{E}\!\left[\left\|Z\right\|_{\mathrm{F}}^{2}\right]-\left\|\mathbb{E}\!\left[Z\right]\right\|_{\mathrm{F}}^{2}.$$

With these definitions, Chebyshev’s inequality (which is just Markov’s inequality applied to the real-valued random variable $\left\|Z-\mathbb{E}\!\left[Z\right]\right\|_{\mathrm{F}}^{2}$) holds for matrix-valued random variables.

We will also need MacDiarmid’s bounded difference inequality [34].

In order to amplify the success probability of probabilistic estimators, we will often use the following version of the “powering lemma” from [25]. Since our version is slightly more general than the original version (our version applies to complex-valued random variables as well), we include a proof for completeness.

We first define query access and (approximate) sampling access to a vector.

We now define the concept of (approximate) sampling-and-query access, which is the main definition of this paper.

Knowledge of $\left\|u\right\|$ is required since this norm cannot be easily computed from $\mathsf{Q}(u)$ and $\mathsf{S}_{\varepsilon}(u)$. The literature on sampling methods for linear algebra, as well as recent works on dequantization [5, 39], have used $\mathsf{SQ}(u)$ (i.e., Definition 3 with $\varepsilon=0$).¹⁴¹⁴14We also mention the work by Cotler, Huang and McClean [8], which compares $\mathsf{SQ}(u)$ with quantum access to a quantum state proportional to $u$ and shows that $\mathsf{SQ}(u)$ can in some cases be more powerful due to the query access allowed in the classical setting (this advantage naturally disappears if we allow quantum algorithms to use $\mathsf{Q}(u)$). We are not aware of any prior work that introduces a similar definition with $\varepsilon>0$.

We define two generalizations of our central definition (Definition 3) that will be useful to state our results.

The first one is a slight generalization of Definition 3, which only requires to get an approximation of $\left\|u\right\|$ and only requires randomized procedures for sampling and norm estimation. This generalization is needed since in some applications we will not be able to implement sampling with probability one or compute exactly the norm of the vector.

The second one is a more significant generalization based on the concept of oversampling from [13]. This is a generalization of the definition from [5], which corresponds to the case $\varepsilon=0$. This generalization is convenient since (as we will show later) it has several useful closure properties.