A parallel quantum eigensolver for quantum machine learning

The past decade has witnessed the great successes of machine learning in many areas [1]. However, with the end of Moore's law and the rapidly increasing demands for machine learning, it is necessary to develop new computing machines [2–4]. Quantum computing is a potential candidate that has exceeded modern supercomputers in the specific tasks of random circuit sampling and boson sampling [5, 6]. Quantum machine learning (QML) implies that quantum computing may boost machine learning tasks using quantum devices [7–11]. In QML, solving eigenproblems is often encountered in practical applications [12].

Many QML algorithms based on solving eigenproblems are designed, such as quantum recommendation system (QRS) [13], quantum linear solver (QLS) [14–18], quantum singular value thresholding (QSVT) [19], quantum support vector machine (QSVM) [9], and quantum principle component analysis (QPCA) [20], but only some of them are realized in experiments, such as QPCA on $4\times 4$ matrices [21, 22] and QSVM on $2\times 2$ training data [23]. QLS models are implemented in superconducting qubits [15, 24, 25], nuclear magnetic resonance (NMR) [26, 27], and photonics [28, 29]. One of the main challenges is the implementation of quantum phase estimation (QPE) because it costs many extra qubits to ensure high precision [30]. Some variational quantum algorithms (VQAs) are also proposed to explore quantum advantages on Noisy Intermediate-Scale Quantum (NISQ) devices, such as variational quantum eigensolver and variational quantum linear solver [15, 31–37]. There are also some works that utilize linear combination of unitaries [38] or block encodings to solve linear algebra problems [39–41]. Recently, quantum resonance transition (QRT) that merely needs two ancillary qubits and achieves squared acceleration over QPE has been proposed to solve the eigenproblem of physical systems [42]. It is further developed to estimate the energy spectrum of the $\textrm{H}_\textrm{2}\textrm{O}$ molecule on NMR [43] and to construct QPCA on the nitrogen-vacancy center [22].

In this work, inspired by QRT, we construct quantum multi-resonant transitions and propose parallel quantum eigensolver (PQE) by designing its Hamiltonian dynamics. PQE directly solves a eigenproblem of the projective summation of a vector b in the partial eigenspace of the hermitian matrix A in the given weights, instead of individually solving each eigenvector. We develop PQE-based collaborative filtering QRS and QLS in the time complexity $\mathcal{O}(\kappa d R\textrm{ploylog}N/\epsilon)$ with the condition number κ, dimension N and sparsity d, the number of non-degenerate eigenvalues R of matrix A, accuracy and some assumption. We implement a movie recommendation experiment on a four-qubit NMR quantum processor, where a $4\times 4$ movie-scoring table (four users and four movies [44]) is considered. The experimental results show that the recommender can predict the scores on unwatched movies and accurately recommend movies to user that he/she might like. Finally, we numerically simulate the applications of PQE in data classification and image completion.

Given a matrix A, an N-dimensional vector b , and a function $f(\lambda_j)$ where λ_j 's are partial or all eigenvalues $\{\lambda_1, \lambda_2, {\ldots}, \lambda_R\}$ of A, how to obtain the solution

$$\begin{equation} \boldsymbol{x} = \sum_{k = 1}^{R}f\left(\lambda_k\right) \left(\boldsymbol{v}_k \cdot \boldsymbol{b}\right) \boldsymbol{v}_k , \end{equation} \tag{ 1 }$$

where R is the number of eigenvalues to be considered and v _k is the eigenvector of A. It is the projective summation of b in the partial eigenspace of A according to the given weight $f(\lambda_k)$. $f(\lambda_k)$ depends on specific problems. For instance, QLS can be developed by $f(\lambda) = \lambda^{-1}$ [45]. Hence, studying how to solve this eigenproblem will generate some QML applications. One direct idea is to use existing quantum eigensolvers [20, 42, 46, 47] to individually solve each λ_k and v _k , and sum directly in equation (1), but it consumes lots of resources and cannot maintain the correct relative phase between individual $|v_k\rangle$ [48, 49]. PQE can directly solve this eigenproblem without the need of diagonalization. It needs one probe qubit, an ancillary register with $r = \lceil\log_2R\rceil$ qubits, and a working system with $n = \lceil\log_2 N\rceil$ qubits. It includes critical steps as follows,

(a) Hamiltonian design. The $(1+r+n)$-qubit Hamiltonian is constructed as

$$\begin{equation} \begin{aligned} \mathcal{H}_{\textrm{PQE}}& = \frac{1}{2}\sigma_z \otimes I_2^{\otimes \left(r+n\right)} + |1\rangle\langle 1|\otimes H_\lambda \otimes I_2^{\otimes n} \\ &\quad + |1\rangle\langle 1|\otimes I_2^{\otimes r} \otimes A + c \sigma_x \otimes S \otimes I_2^{\otimes n}. \end{aligned} \end{equation} \tag{ 2 }$$

c is the coupling parameter, I₂ is a $2\times 2$ identity and $\sigma_{x, y, z}$ are the Pauli operators. The four parts of equation (2) describe the Hamiltonians of probe qubit, ancillary register, working system, and the interaction between them, respectively. H_λ and S are diagonal matrices. H_λ can trigger multi-resonant transitions by encoding the eigenvalues of resonant energy levels into the diagonal elements $h_k = 1-\lambda_k$. The element s_k in S controls the strength of the kth transition. Figure 1 illustrates the mechanism design.

Zoom In Zoom Out Reset image size

(b) Dynamical evolution. The system starts from $|\Psi\rangle_\mathrm{in} = |0\rangle \sum_{k = 0}^{R-1}|k\rangle | b\rangle/\sqrt{R}$ and evolves under $\mathcal{H}_{\textrm{PQE}}$, where $| b\rangle = \boldsymbol{b}/||\boldsymbol{b}||$ and $|k\rangle$ is the computational basis of the ancillary register. As shown in figure 1. $\mathcal{H}_{\textrm{PQE}}$ is a sparse block matrix. In the subspace of $\{|0\rangle|k\rangle|v_k\rangle, |1\rangle|k\rangle|v_k\rangle\}$, $H_{00} = H_{11} = 1/2$ and $H_{01} = H^{*}_{10} = cs_k$. Its evolution corresponds to a rotation from $|0\rangle|k\rangle|v_k\rangle$ to $|1\rangle|k\rangle|v_k\rangle$ with the strength of cs_k . There is still a small transition from $|0\rangle |k\rangle |v_j\rangle$ to $|1\rangle |k\rangle |v_j\rangle$ for $|v_j\rangle (j \neq k)$. Ignoring this effect for the moment, $|\Psi\rangle_\mathrm{in}$ evolves into $|\Psi_1(t)\rangle$,

$$\begin{equation} \begin{aligned} \frac{e^{-it/2}}{i\sqrt{R}} |1\rangle \sum_{k = 0}^{R-1}\sin\left(cs_k t\right)\langle v_k| b\rangle |k\rangle |v_k \rangle +\alpha_{\overline{1}} |\Psi\rangle_\mathrm{else}. \end{aligned} \end{equation} \tag{ 3 }$$

After Hadamard gates are applied on the ancillary register, the final state $|\Psi_2(t)\rangle$ is

$$\begin{equation} \begin{aligned} \frac{e^{-it/2}}{iR} |1\rangle|0\rangle ^{\otimes r} \sum_{k = 0}^{R-1}\sin\left(cs_k t\right)\langle v_k| b\rangle |v_k \rangle +\alpha^{^{\prime}}_{\overline{10}} |\Psi^{^{\prime}}\rangle_\mathrm{else}. \end{aligned} \end{equation} \tag{ 4 }$$

$\alpha_{\overline{1}}$ and $\alpha^{^{\prime}}_{\overline{10}}$ are the factors in uninterested subspaces.

(c) Measurements. We set $t = \frac{1}{c}$ and $s_k = \arcsin\frac {f(\lambda_k)}{f_{\max}}$ with $f_{\max, \min} = \max,\min_{\lambda_k} |f(\lambda_k)|$. Then we measure the work system in the ancillary subspace of $|1\rangle|0\rangle^{\otimes r}$. The final result is

$$\begin{equation} \boldsymbol{x}^{^{\prime}}_{\textrm{PQE}} = \frac{e^{-i/2c}}{iRf_{\max}} \sum_{k = 0}^{R-1} f\left(\lambda_k\right)\langle v_k| b\rangle |v_k \rangle. \end{equation} \tag{ 5 }$$

It differs from equation (1) by a constant factor. Compared with quantum single-resonant transitions [42, 50, 51], we add an r-qubit register that encodes the eigenvalues, and we design the Hamiltonian $\mathcal{H}_{\textrm{PQE}}$ and its dynamics in PQE, such that PQE can trigger multi-resonant transitions simultaneously. Thus, it can directly obtain equation (5) without diagonalizing the matrix A. Here, the main steps in the complexity analysis are briefly presented. First, the success probability of the measurement at the end of the quantum circuit can be easily calculated. It is,

$$\begin{align*} P_m = \frac{1}{R^2f_{\max}^2} \sum_{k = 0}^{R-1}\left[f\left(\lambda_k\right)b_k\right]^2 \unicode{x2A7E} \frac{1}{R^2}\frac{f_{\min}^2}{f_{\max}^2} \sum_{k = 0}^{R-1}b_k^2 \end{align*}$$

where $b_k = \langle v_k| b\rangle$. This probability can be increased to a constant by $\mathcal{O}(1/\sqrt{P_m})$ times amplitude amplification [52]. Second, the query complexity of implementing the evolution of a d-sparse Hamiltonian H is $\mathcal{O}(td||H||_{\textrm{max}})$ pre amplitude amplification. Third, if the error is defined as the distance between the PQE result $\boldsymbol{x}^{^{\prime}}_{\textrm{PQE}}$ and the ideal result x , the magnitude of error introduced by the off-resonant effect is $\epsilon \sim c/|\lambda_j-\lambda_k|$. When the interval between the eigenvalues is $\Omega (1 / \textrm{polylog}N)$, the total time complexity of PQE will be $\mathcal{O}(||H||_{\textrm{max}}\kappa d R\beta^{-0.5}\textrm{polylog}N/\epsilon)$ with $\kappa = f_{\max}/f_{\min}$ and $\beta = \sum_{k = 0}^{R-1}b_k^2$. Similar to the previous [14], we also assume that $||H||_{\textrm{max}}$ is bounded by some constant. We stress that the complexity of PQE depends on the specific problems, for example, $\beta\sim 1$ and the final complexity will be $\mathcal{O}(\kappa d R\textrm{polylog}N/\epsilon)$ in the QRS and QLS. More details can be found in supplement material [45]. Next, we demonstrate its application in QRS [53].

In the various recommendation systems, collaborative filtering is the mainstream model which uses the similarity between users or items to recommend the target items to the users [54]. The item-based collaborative filtering recommendation system, as an example, contains three steps shown in figure 2. First, collect the scoring table B with M users and N items. Second, compute the similarity between the ith and jth items by dot product $\textrm{sim}(i, j) = B_{:i}\cdot B_{:j}$ with column-normalized B. Here, singular value decomposition [52, 55] is adopted to decompose B into a low-dimensional space and then compute the dot product matrix by $\textrm{sim} = \sum_{k = 1}^R\Lambda_k^2 V_{:k} V_{:k}^T$ with the first R maximum singular values $\Lambda_k$. V is the right singular matrix of B. Third, calculate the prediction of a new-user vector b on the jth item by

$$\begin{equation} \begin{aligned} W_j = \sum_{i = 1}^{N} b_i \cdot \textrm{sim}\left(i,j\right) = \sum_{k = 1}^R\Lambda_k^2 \left(V_{:k}\cdot \boldsymbol{b}\right)V_{jk}. \end{aligned} \end{equation} \tag{ 6 }$$

The system will recommend the item with the highest score to the user. Similar procedure can be made for the user-based collaborative filtering recommendation system.

Zoom In Zoom Out Reset image size

Obviously, equation (6) is the solution of equation (1), so PQE can realize the item-based QRS by setting $A = B^{\dagger}B$ and $f(\lambda_k) = \lambda_k$. The unitary operator $e^{-iB^\dagger Bt}$ can be implemented by the method in [9]. Besides, we can also set

$$\begin{align*} A = \left[\begin{array}{cc} 0 & B \\ B^\dagger & 0 \end{array}\right] \end{align*}$$

to construct the Hermitian matrix containing the information of B in QRS with input state $|b\rangle\propto[\textit{b},\ \boldsymbol{0}]^\dagger$. Figure 2 presents the workflow schematic for the PQE-based QRS by taking the item-based case as an example. $\mathcal{H}_{\textrm{PQE}}$ can be directly performed if the eigenvalues of A are known. Otherwise, some quantum eigensolvers are needed to obtain the first R larger-weight eigenvalues before that. Here, we propose an improved QRT to finish this task [45].

A normalized state vector $|b\rangle = \sum_i b_i|i\rangle/||\boldsymbol{b}||$ encodes a new-user vector b . We feed $\mathcal{H}_{\textrm{PQE}}$ and $|b\rangle$ into PQE and execute it. In the ancillary subspace of $|1\rangle|0\rangle ^{\otimes r}$, measuring the work system will output the result equation (6) by a constant factor. Considering that scores are positive values and only the probabilities of wave functions can be measured, we use an experiment-friendly definition as our result by $W^{\textrm{PQE}}\sim |\boldsymbol{x}^{^{\prime}}_{\textrm{PQE}}| = \sqrt{\textrm{Diag}(\rho_f^{|10\rangle})}$. $\textrm{Diag}(\rho_f^{|10\rangle})$ are the diagonal elements of final density matrix ρ_f in the ancillary space of $|1\rangle|0\rangle ^{\otimes r}$. For ensemble quantum devices, like NMR $\textrm{Diag}(\rho_f^{|10\rangle})$ can be inferred by measuring the operators $\{I, \sigma_z\}^{\otimes n}$ [56]. For quantum devices with projective measurements such as superconducting circuits [57], W^PQE is easily obtained on the computational basis of work qubits. Finally, we can recommend the item with the highest score to the new user. The user-based case can be similarly obtained using $A = BB^{\dagger}$ and row-normalized B.

As a verification, we implement a movie recommendation experiment on a four-qubit NMR Here, we consider a small-scale user-movie scoring table with four users and four movies that is selected from GroupLens [44],

$$\begin{align} B = \left(\begin{array}{ccccc} &\textrm{M-1} & \textrm{M-2} & \textrm{M-3} & \textrm{M-4} \\ \textrm{U-1} & 5 & 4 & 2 & \\ \textrm{U-2} &5 & 5 & & \\ \textrm{U-3} &1 & & 4 & 4 \\ \textrm{U-4} & & & 5 & 5\end{array}\right).\nonumber \end{align}$$

The four-qubit sample is unlabeled iodotrifluroethylene (C₂F₃I) molecule dissolved in d-chloroform, in which the half-spin nucleus (one ¹³C and three ¹⁹F) in the sample are used as the qubits. The structure and properties of the molecule are shown in figure 3(a), where one ¹³C and three ¹⁹F form a four-qubit all-to-all connected quantum processor. ¹³C has both the abilities of individual address and readout, while ¹⁹F nucleus lack the individual readout. This problem can be overcome by transferring ¹⁹F channel to ¹³C channel using SWAP operation in experiments. Under weak coupling approximation, the internal Hamiltonian of this 4-qubit system is

$$\begin{equation} \mathcal{H}_\mathrm{int} = -\sum\limits_{j = 1}^4 {\pi \nu _j } \sigma_z^j + \sum\limits_{j < k, = 1}^4 {\frac{\pi}{2}} J_{jk} \sigma_z^j \sigma_z^k, \end{equation} \tag{ 7 }$$

where ν_j is the chemical shift and J_jk is the J-coupling strength between jth and kth nuclei. The dynamic of the spins is controlled by designing the shaped radio-frequency pulses, which can be used to realize the evolution of quantum multi-resonant Hamiltonians.

Zoom In Zoom Out Reset image size

First, we implement QRT to search for the first two eigenvalues that are necessary for the subsequent procedure. It needs one probe qubit (¹³C) and two work qubits (first two ¹⁹F). The Hamiltonian is $\mathcal{H}_{\textrm{QRT}} = \omega_x |0\rangle\langle 0|^{\otimes 3}+ |1\rangle\langle 1| \otimes A + c \sigma_x\otimes I_2^{\otimes 2}$. Starting from $|0\rangle^{\otimes 3}$, we first roughly locate the resonant positions by setting c = 0.1 and sweeping ω_x from 1 to 3 with the step $\bigtriangleup \omega_x = 0.1$ and then accurately find the eigenvalues by fixing c = 0.01 and sweeping ω_x near resonance peaks. We optimize the shaped pulses by the gradient ascent pulse engineering technique [58, 59] and implement $U_{\textrm{QRT}} = e^{-i\mathcal{H}_{\textrm{QRT}}\tau}$ with $\tau = 2/c$ for each pair of c and ω_x . The width of the shaped pulses is 20 ms. The probability $P_{|1\rangle}$ of probe qubit is obtained by $P_{|1\rangle} = 1-\textrm{Tr}(\rho_f\sigma_z^1)/8$ [45]. Figure 4 presents the experimental results of sweeping ω_x for $A = B^{\dagger}B$ (item-based) and $A = BB^{\dagger}$ (user-based). For $A = B^{\dagger}B$, there are two resonance peaks corresponding to two eigenvalues $\lambda_1^\mathrm{exp} = 2.263~(\lambda_1^\mathrm{th} = 2.264)$ and $\lambda_2^\mathrm{exp} = 1.700~(\lambda_2^\mathrm{th} = 1.699)$. $\lambda_1^\mathrm{exp} = 2.316~(\lambda_1^\mathrm{th} = 2.316)$ and $\lambda_2^\mathrm{exp} = 1.650~(\lambda_2^\mathrm{th} = 1.649)$ are estimated for $A = BB^{\dagger}$. The remaining smaller eigenvalues are outside the frequency range of the scan.

Zoom In Zoom Out Reset image size

Second, we perform PQE with $\lambda_1^\mathrm{exp}$ and $\lambda_2^\mathrm{exp}$. It needs one probe qubit (¹³C), one register qubit (¹⁹F₁) and two work qubits (the remaining ¹⁹F). To prepare the initial state $|b\rangle$, it usually needs some quantum devices for a general $|b\rangle$ like quantum random access memory. Fortunately, the initial state $|b\rangle$ is usually sparse enough in the QRS to be prepared from a product state through a few gates. In the item-based recommendation, we consider the recommendation task for two new users U-5 and U-6 with vectors $|b\rangle^{\textrm{U-5}}\propto[5,0,1,0]^T$ and $|b\rangle^{\textrm{U-6}}\propto[1,0,5,0]^T$ respectively. They have watched M-1 and M-3 but did not watch M-2 and M-4. The task is to predict their scores on M-2 and M-4 and recommend the movie to them. In the user-based case, we recommend one of two new movies M-5 and M-6 to the users. Their scoring vectors are $|b\rangle^{\textrm{M-5}}\propto[5,0,2,0]^T$ and $|b\rangle^{\textrm{M-6}}\propto[2,0,5,0]^T$. U-1 and U-3 watched and scored new movies, but U-2 and U-4 have not watched them. The task is to predict the scores of U-2 and U-4 on new movies and recommend the movie to them.

As shown in figure 3(b), the experiment includes three procedures. (i) Preparing the initial state. Single-qubit rotations are applied on ¹⁹F₁ and ¹⁹F₂ to prepare the state $|0\rangle|+\rangle|b\rangle$ from $|0\rangle^{\otimes 4}$. The shaped pulses are used realize them in experiments. (ii) Realizing the evolution operators. The shaped pulse is optimized to implement the evolution $e^{-i\mathcal{H}_{\textrm{PQE}}/c}$, similar to the above experiments. The last Hadamard gate on ¹⁹F₁ is decomposed into single-qubit rotations. Totally, we engineer the implementation of the quantum circuit with a 15 ms shaped pulse. (iii) Obtaining recommendation results. By using $\pi/2$ readout pulses and collecting the NMR signal spectrum, we measure the expectation values of $\{I, \sigma_z\}^{\otimes 4}$ to infer the diagonal elements $\textrm{Diag}(\rho_f^{|10\rangle})$ [45]. Before the measurement, we also perform the amplitude amplification to increase the probability $P_{|1\rangle|0\rangle}$ [45, 60]. It is also realized with a 20 ms shaped pulse. Finally, the movie with the highest score will be recommended to the user. Equation (6) shows that W will give information about similarities between movies or users when the input $|b\rangle$ is a simple product state. Here, we also individually feed product states $|b\rangle = |00\rangle$ and $|b\rangle = |10\rangle$ into both item-based and user-based QRS experiments. Their results $W^{\textrm{PQE}}_{|00\rangle}$ and $W^{\textrm{PQE}}_{|10\rangle}$ respectively provide the similarities of M-1 and M-3 with the remaining movies (item-based) and the similarities of U-1 and U-3 with the remaining users (user-based).

Figure 5 presents the experimental results of both item-based and user-based QRS. The results agree well with the numerical simulations, verifying the feasibility of the PQE. In the item-based QRS, the user U-5 has watched the movies M-1 and M-3 and rated the former with a higher score. For unwatched movies, the recommendation system rated M-2 higher than M-4. Thus, it will recommend M-2 to U-5. It is an accurate recommendation because both M-1 and M-2 are action movies [44]. Therefore U-5 is likely to appreciate M-2. Experimental results also support this recommendation, of which M-2 has a higher similarity with M-1 than M-4. Based on a similar analysis, it will recommend M-4 to U-6. In the user-based QRS, U-1 and U-3 have watched M-5, but U-2 and U-4 did not watch it. The results show that U-2 will rate M-5 a higher score than U-4 after watching M-5. Thus M-5 will be recommended to U-2. It makes sense that U-1 rates M-5 highly because the taste between U-1 and U-2 is similar. Namely, they may both enjoy action movies. With the same logic, U-4 might prefer M-6 since U-3 rated M-6 highly.

Zoom In Zoom Out Reset image size

The deviation between experimental and ideal diagonal elements is 1.56%, mainly from the imperfections of pulses and decoherence. The simulation shows that their contribution is approximately 1.49% [45]. The accuracy of the found eigenvalues may affect the performance of PQE. In QRT experiments, the found eigenvalues have a deviation of 10⁻³ from the real one, and $c = 10^{-2}$ in PQE. In such a setting, the final state of PQE has a fidelity of over 0.998 [45].

Our work can be extended to solve other problems in machine learning, like classification and image completion, by setting different functions f in PQE. We simulate the application in data classification by using our method as a quantum linear solver and image completion by using our method to realize singular value threshold (SVT).

For data classification, given the training data with M points $\{(\vec{z}_j, y_j)| \vec{z}_j \in \mathbb{R}^N, y_j = \pm1, j = 1, {\ldots}, M\}$, the support vector machine (SVM) classifies a new point $\vec{z}_\textrm{new}$ into one of two classes. Here, $y_j = \pm1$ is the label of the points. The critical step of predicting sample classification is to solve the linear equation,

$$\begin{equation} \left(\begin{array}{l} \eta_0 \\ \vec{\eta} \end{array}\right) = A^{-1}\left(\begin{array}{l} 0 \\ \vec{y} \end{array}\right), A = \left(\begin{array}{cc} 0 & \overrightarrow{1}^{T} \\ \overrightarrow{1} & \mathcal{K}+\gamma^{-1} I \end{array}\right). \end{equation} \tag{ 8 }$$

$\mathcal{K}$ is the kernel matrix with the element $\mathcal{K}_{jk} = \vec{z}^T_j\cdot \vec{z}_k$. $\vec{y} = \left(y_{1}, \ldots, y_{M}\right)^{T}$ and $\overrightarrow{1} = (1, \ldots, 1)^{T}$. PQE solves the above equation by constructing the Hamiltonian of the composite system with one probe qubit, an $r = \lceil \log_2 R\rceil$-qubit register, and an n-qubit work system. The Hamiltonian is designed as equation (2) with the function form $f(\lambda) = \lambda^{-1}$. The system starts from the initial state $|\Psi\rangle_{\textrm{in}} = |0\rangle\otimes\sum_{k = 0}^{R-1}|k\rangle\otimes |b\rangle$. Here, $|b\rangle$ is the normalized vector of $(0, \vec{y})$. The Hadamard gates are applied on the register after time dynamics with $t = 1/c$. Measuring the state of work qubits in the ancillary subspace of $|1\rangle|0\rangle^{\otimes r}$ gives the result x_work. The solution of equation (8) can be written as,

$$\begin{equation} \left(\begin{array}{l} \eta_0 \\ \vec{\eta} \end{array}\right)_{\textrm{PQE}} = \mathbb{N} \cdot \frac{iRf_{\max}}{e^{-i/2c}} \cdot x_{\textrm{work}}. \end{equation} \tag{ 9 }$$

Here, $\mathbb{N}$ is the normalized factor of $(0, \vec{y})^T$. This solution can be further used to classify new data into one of two classes [45]. Some classification results are placed in figure 6.

Zoom In Zoom Out Reset image size

For image completion, SVT is one of the common methods in matrix completion. For a sampled matrix $B\in \mathbb{R}^{N \times M}$, it has the singular value decomposition $B = \sum_{k = 1}^{R}\Lambda_k\vec{u}_k \vec{v}_k^\dagger$, where $\vec{u}_k$ and $\vec{v}_k$ are left and right singular vectors corresponding to singular value $\Lambda_k$. The task of SVT is to recover a new matrix B' with the eigenvalue threshold $\Lambda_{\tau}$,

$$\begin{equation} \begin{aligned} B^{^{\prime}} = \sum_{k = 1}^{R} f_\tau\left(\Lambda_k\right) \vec{u}_k \vec{v}_k^\dagger, ~ f_\tau\left(\Lambda\right) = \max\left(\Lambda-\Lambda_{\tau},0\right). \end{aligned} \end{equation} \tag{ 10 }$$

To solve this problem, we need to design a Hamiltonian, which is a little different from equation (2),

$$\begin{equation} \begin{aligned} \mathcal{H}_{\textrm{SVT}} = &\left(\frac{1}{2}\sigma_z \otimes I_2^{\otimes r+n_1} + |1\rangle\langle 1|\otimes H_\Lambda\otimes I_2^{\otimes n_1} \right. \\\nonumber & +\left. |1\rangle\langle 1|\otimes I_2^{\otimes r} \otimes A + c \sigma_x \otimes S \otimes I_2^{\otimes n_1}\right)\otimes I_2^{\otimes n_2}. \end{aligned} \end{equation}$$

Here, $A = BB^\dagger$. S and $H_\Lambda$ are the diagonal matrices with the elements $s_k = \arcsin\frac{f_\tau(\Lambda_k)}{\Lambda_k}$ and $h_k = 1-\Lambda_k^2$, respectively. The other next steps are the same as QRS. The initial state denotes the input matrix B, and the completed matrix B' will be encoded in the final state. Figure 6 presents the simulated results, and more details can be found in supplemental information [45].

In this section, we first make some discussions related to QRT. The time complexity of QRT is $\mathcal{O}(R d\textrm{log}N /\delta)$ with the precision δ of eigenvalues [45]. It is worth noting that QRT is not the only option to solve the eigenvalues. There are some quantum approaches here that can solve the eigenvalues, such as near-optimal ground state preparation using a method called block coding [61], but this approach seems to solve directly only for the ground state and requires more ancilla qubits. In this work, we propose a modified QRT and use it to find eigenvalues, as it requires only one ancilla qubit and can easily obtain several large eigenvalues, which is very friendly for our experiments.

We second discuss the complexity and error source of the PQE. There is still a small transition from $|0\rangle |k\rangle |v_j\rangle$ to $|1\rangle |k\rangle |v_j\rangle$ even if $c\ll|\lambda_j - \lambda_k|$, which introduces an error and this error $\epsilon = ||x_{\textrm{PQE}}-x||$ is related with $c/|\lambda_j - \lambda_k|$. When the interval between the eigenvalues is $\mathcal{O} (1 /\textrm{ploylog}N)$, the total complexity of PQE is $\mathcal{O}(\kappa d R\beta^{-0.5}\textrm{polylog}N/\epsilon)$ with sparsity d, $\beta = \sum_{j = 1}^{R}|b_j|^2$, and $\kappa = f_{\textrm{max}}/f_{\textrm{min}}$ [45]. Here, we assume that $\|H\|_{\max}$ is bounded by some constant. Hence, the scaling behavior of PQE depends on the problem studied. PQE-based QRS and QLS has the complexity $\mathcal{O}(\kappa d R\textrm{polylog}N/\epsilon)$. Meanwhile, PQE still works when the degenerate eigenvalues exist because resonant transitions are triggered according to eigenvalues instead of eigenstates [42].

We third make a comparison with the previous QRS and QLS. Different from previous QRS [13], our work is based on a collaborative filtering recommendation system rather than low-rank approximate matrix completion. Besides, it solves the eigenproblem using the PQE proposed in this paper rather than quantum phase estimation, which uses fewer qubits than the previous work. To our knowledge, the best collaborative filtering QRS before PQE has linear scaling in N [62]. Thus, PQE for collaborative filtering QRS can reduce the complexity to polylogarithmic scaling in N when $R = \mathcal{O}(\textrm{log}N)$. Compared with the HHL method whose complexity is $\mathcal{O}(\kappa^2 d^2 \textrm{log}N/\epsilon)$, PQE has linear scaling in κ and d, and it may achieve polynomial speedup when $R = \mathcal{O}(\textrm{log}N)$. For QLS, it can be solved by some works based on linear combinations of unitaries [18, 63]. They need to calculate $\mathcal{O}(\log (\kappa/\epsilon)\kappa/\epsilon \times \kappa \log (\kappa/\epsilon))$ different overlaps with $\mathcal{O}(\kappa \log (\kappa/\epsilon))$ evolution time, or compute a linear combination of so many unitary operators. In contrast to them, the number of ancilla qubits $1+\lceil\log_2 R\rceil$ is independent of κ and ε in PQE. At the same time, PQE does not require the measurement of many overlaps or summing over a lot of terms.

Last, we also compare the PQE with VQAs-based NISQ methods. VQAs have different structures and complexities for different problems and currently suffer from some challenges in trainability and efficiency [32–34]. PQE has a definite structure and complexity, and it does not need optimization, but its essential constituent is the simulation of Hamiltonians. There are some efficient methods for realizing Hamiltonian dynamics, such as Taylor-series methods [64–67], product-formula-based methods [68–72], VQA-based methods [73–77], and quantum signal processing [78, 79].

This letter focuses on solving eigenproblems in QML. We propose PQE that can solve the eigenproblems in parallel and further demonstrate the feasibility of PQE by constructing QRS and implementing the movie recommendation experiments. The first experimental realization of QRS benefits from the implementable dynamics and fewer qubit resource requirement. Our framework is friendly to quantum devices and it will pave the way towards more practical applications in QML. Here, we numerically simulate the applications in data classification (PQE is used as QLS) and image completion ( PQE is used to realize QSVT) [45]. The images are accurately recovered even under the low sampling rate of 50%. Moreover, it is expected to develop quantum fidelity estimators without full state tomography [80, 81] and the SWAP test [82], because the overlap $\langle v_k | b\rangle$ is computable in PQE. It is also interesting to explore quantum chemistry simulations [83].

This work is supported by Guangdong Basic and Applied Basic Research Foundation (2022B1515020074), the National Natural Science Foundation of China (12275117), Guangdong Provincial Key Laboratory (2019B121203002), Beijing Nova Program 20230484345, Shenzhen Science and Technology Program (KQTD20200820113010023), and Tsinghua University Initiative Scientific Research Program, Beijing Advanced Innovation Center for Future Chip (ICFC).

All data that support the findings of this study are included within the article (and any supplementary files).