Quantum Quantitative Trading: High-Frequency Statistical Arbitrage Algorithm

Quantitative trading is an integral part of financial markets with high calculation speed requirements, while no quantum algorithms have been introduced into this field yet. We propose quantum algorithms for high-frequency statistical arbitrage trading in this work by utilizing variable time condition number estimation and quantum linear regression.The algorithm complexity has been reduced from the classical benchmark O(N^2d) to O(sqrt(d)(kappa)^2(log(1/epsilon))^2 )). It shows quantum advantage, where N is the length of trading data, and d is the number of stocks, kappa is the condition number and epsilon is the desired precision. Moreover, two tool algorithms for condition number estimation and cointegration test are developed.


I. INTRODUCTION
With the rapid development of quantum computing [1][2][3], the qubits on the chips are up to 53 currently [3], and it will extend beyond 100 soon in the roadmap of quantum systems based on superconductivity.Hence, quantum computing shows the potential to solving practical problems, such as chemistry [4][5][6], materials [7], drug design [8], and et al.
Quantitative trading is an essential field of finance, and statistical arbitrage is a mainstream approach of quantitative trading taken by most hedge funds [24,25].While lots of classical algorithms for quantitative trading have been proposed [26][27][28][29], and traditional hardware techniques including infrared communication and Field Programmable Gate Array have been employed over the years [30,31], still the requirement for speed cannot be satisfied when implementing those complicated statistical methods, especially in the quicker-take-all situation of high-frequency trading(HFT) whose need of computing speed is crucial [32].In statistical arbitrage, one needs to find a potential cointegrated pair via many linear regressions and cointegration tests involving a huge matrix of historical data.For example, in U.S. stock markets, the problem size can exceed N = 10 7 and the complexity is 10 15 (see section VI for details) which is very hard to calculate by classical computers.For this problem, quantum computation might provide an effective solution.
In this article, quantum algorithms applied to statistical arbitrage strategy are proposed.It consists of two subroutines: the first one is the Variable Time Preselection Algorithm(VTPA) that will help to find, with high probability, the potential comovement out of securities and portfolios.The second one is the Quantum Cointegration Test Algorithm(QCTA) that focuses on the efficient verification of cointegrated pairs, which is quite valuable in statistics but has not been achieved via quantum computation ever before.The classical benchmark to achieve the preselection procedure is by matrix factorization with complexity O(N 3 ) [33], while our algorithm's complexity is O( )) where d is the number of stocks usually much less than time length N and κ 0 is the condition number.Moreover, an efficient tool named Quantum Condition Number Comparison Algorithm (QCNCA) used to probe a matrix's condition number is proposed, and it can be applied to many other domains.
The structure of this article is as follows: After giving the preliminaries in section II, the global structure and main results of our work are shown in section III.The details of VTPA and QCT are described in sections IV and V, respectively, followed by a discussion on complexity and quantum advantage in section VI.

II. PRELIMINARIES
Since different domains, including quantum computing, statistics, and finance, are covered while the readers may not be familiar with one of them, related preliminaries are introduced in detail.

A. multicollinearity
In this subsection, a brief introduction about multicollinearity and condition number will be given, helping to understand the first preselection algorithm.In statistics, multicollinearity refers to a situation in which some of the explanatory variables in a multiple regression model are highly linearly related.
In numerical analysis, to detect and measure the seriousness of the multicollinearity problem, the condition number κ is introduced.Given problem f , it is generally a measurement to describe the change of the output value divided by the change of the input variable x: In the case of matrices, the condition number associated with the linear equation Ax = b release the dependence of accuracy on the input data.Specifically, the condition number of normal matrix A is It should be emphasized that condition number is a property of the matrix itself and does not depend on the algorithm or accuracy of the computer used.Hence, both classical and quantum computers have a common problem to solve an ill-conditioned (high condition number) linear system of equations.The larger the condition number, the more ill-conditioned the matrix is, and the algorithm complexity will increase very quickly.

B. Cointegration
In this subsection, some statistical concepts and facts about stochastic process and time series analysis are provided.Following those, an explicit demonstration is also given on the relationship between multicollinearity and cointegration, which may be confuse some readers.
A (weakly) stationary time series, x t , is a finite variance process with an unconditional joint probability distribution.Thus it does not change when shifted in time: (i) the mean value function µ xt = E(x t ) is constant and (ii) the covariance function γ x (s, t) = E[(x s −µ s )(x t −µ t )] depends on s and t only through their difference|s − t|.In autoregressive-moving average models of unknown order, to test whether a given time series denoted as Y t is stationary or not, the Augmented Dickey-Fuller(ADF) unit root test may be employed [34].
Cointegration (multi-cointegration) is a relevant statistical property of two or more time series which are individually integrated of order d while their combination is integrated of order less than d.Here the order of integration is a summary statistic denoting the minimum number of differences taken to obtain a covariancestationary series.Without loss of generality, d = 1 is assumed in this article.Under different financial hypotheses, there are mainly three kinds of cointegration tests: the Engle-Granger test [35], the Johansen test [36], and the Phillips-Ouliaris test [37].In our work, Engle-Granger two-step method is used as the most popular and famous one: Suppose that x i t are non-stationary and integrated of order d=1, then a linear combination ût = β i x i t is expected to be stationary for some specific coefficient of β i .In the general case that β i is not decided yet, some estimation must be made first, usually by ordinary least squares regression.Next, the stationarity test will be implemented on the residuals ût .It is a regression on ût , and the lagged residuals ût−1 are included as a regressor: Here α and β are the intercept and the coefficient on the time trend, respectively, and p denotes the lag order of the autoregressive process to be decided.The unit root test is carried out under the null hypothesis γ = 0.The test statistic to be computed is .
What follows is a comparison with the Dickey-Fuller distribution critical value table [38].Whenever such a cointegrated stock portfolio is found, the linear combination is expected to have the property of mean-reverting and use the statistical arbitrage.

C. Quantum Linear Regression
Quantum linear regression is the primary tool of QCT and is introdecd as follows.Wiebe, Braun, and Lloyd (WBL)firstly introduced an algorithm for quantum data fitting [39].Building on Harrow, Hassidim, and Lloyd's (HHL) quantum algorithm for linear systems of equations [40], WBL developed a least-squares estimation using Moore-Penrose Pseudo inverse.WBL's algorithms are mainly suited for data sets whose design matrices are sparse and well-conditioned.Given an N dimension s sparse data matrix, the time complexity is O(log N s 3 κ 6 ǫ −1 ), where the condition number given is κ and the accuracy desired is ǫ −1 .With the technique of quantum principal component analysis(qPCA) and singular value decomposition(SVD) [41], Schuld, Sinayskiy, and Petruccione(SSP) came with an algorithm for prediction based on a linear regression model with least-squares optimization [42].The sparseness condition is removed, and the existence of a low-rank approximation is supposed instead.The time complexity is O(log N κ 2 ǫ −3 ), where an improvement of factor κ 4 is made on the condition number at the cost of worse dependence on accuracy by a factor ǫ −2 .Recently, Guoming Wang presents a new quantum algorithm for fitting a linear regression model using least-squares approach [43].This algorithm builds on Low and Chuang's method for Hamiltonian simulation based on qubitization and quantum signal processing [44,45].Childs, Kothari, and Somma (CKS)'s approach is introduced to inverse the matrix derived from SVD [46].Imposing restrictions on the number of adjustable parameters d, and hence the rank of the design matrix, the gate complexity is O( d 1.5 κ 3 ǫ 2 poly[log 2 ( κ ǫδ )]) with the succeeding probability is at least 1 − ǫ.

III. QUANTUM STATISTICAL ARBITRAGE
Pioneered by Gerry Bamberger [47], statistical arbitrage has developed a lot, and the crux and core are to model the comovement.Following the framework first introduced by Vidyamurthy [27], statistical arbitrage is divided mainly into three key steps: Firstly, two or more securities moved together historically in a formation period should be preselected; secondly, some version of the Engle-Granger cointegration test [35] is taken for verification; thirdly, the spread between them in a subsequent trading period is monitored by some optimal entry/exit thresholds.Since the spread of stocks will revert to its historical mean and, the profit can be made from other traders' irrational behavior by longing the oversold securities and shorting the overbought ones at the same time [26].
In this section, two algorithms solving the quantum statistical arbitrage problem are proposed.One is for the case of fixed condition number threshold; the other is for a fixed number of remained portfolios.The formal statement of the quantum statistical arbitrage problem is as follows: Given historical data of many stocks for a long time interval, our target is to select those stocks that are cointegrated.The algorithm mainly contains two steps: preselect multicollinear stock portfolios from the pool by applying V T P A(p, κ) where V P T A is True if the given portfolio p's condition number is larger than the threshold κ; and then verify whether the preselected portfolio p is cointegrated by implementing QCT (p) to output cointegration flag f and corresponding coefficients β.
Suppose that P = {p} is the portfolio pool, and (p (j) t ) J×T is a portfolio of stocks' historical quote data.Here p (j) t is an element of p as the j th stock's price at time t.The matrix p is of full rank since no perfect linear relation exists in noisy financial market data.The two quantum statistical algorithms work in the standard oracle model, and the matrix is stored in a quantum random access memory(qRAM) [48][49][50].A procedure P x is assumed to perform the map for any j ∈ [1, 2, ..., d] and t ∈ [1, 2, ..., N ], and the price is stored as a bit string in the third register.
In order to derive the desired real symmetric matrix, the strategy of HHL [39,40] is adopted as: x T 0 .Moreover, the norm of the matrix is assumed to satisfy A = 1 without loss of generality since otherwise let A = A A .
If an efficient κ 0 derived from historical data is taken as filter threshold, the following Algorithm 1 is given: As for the case of unknown κ 0 , an even more efficient Algorithm 2 is provided.The basic idea is as follows: since our single-step preselection sub-algorithm can be used for any given κ, a progressive κ preselection procedure can be implemented.Portfolio matrices with small κ will be directly obsoleted in the first several steps until the number of matrices left is small enough, and until then, the quantum cointegration test will be implemented.

Input:
k: portfolio number threshold T : the length of time interval J: the total number of stocks d: the number of stocks in one portfolio P : the portfolio pool p (j) t : the j th stock's price at time t.Output: (p, β) Cointegrated portfolios and cointegration coefficients.
Data Loading. Step Both of the above two algorithms are for statistical arbitrage, and the selection depends on the specific market: if the κ-threshold is stationary, the first algorithm is chosen; otherwise, the second one is preferred.Since the two subroutines are complicated and tool sub-algorithms are developed, they will be introduced in section IV and section V, respectively.

IV. VARIABLE TIME PRESELECTION
In this section, we will explain the main idea of the first part of our work as a variable time quantum algorithm to preselect the stocks that are multicollinear and thus may be cointegrated as needed.
Although ill-conditioned matrices are commonly considered a terrible problem that one should try to avoid, we develop the heuristic idea to detect multicollinearity by searching matrices with small eigenvalues and large condition numbers.QCNCA is developed to determine whether the condition number κ of a given matrix is larger than the threshold κ 0 in subsection A.
Since QCNCA's dependence on κ is quadratic, the technique of variable time quantum algorithm is introduced to accelerate the implementation of matrices selection [51], and then the VTPA is as follows: Theorem 1 Supposing that many different linear systems are given with unknown condition number κ and P j denote the probability that condition number satisfies κ j−1 = 2 j−1 ≤ κ ≤ κ j = 2 j .Then there is an efficient quantum algorithm to preselect matrices with condition numbers κ ≥ κ 0 .The average query complexity is O( As for a uniform probability distribution, the query complexity is O( √ dκ 2 0 log (1/ǫ) 2 ) to determine whether the condition number is larger than κ 0 .
The proofs of correctness and complexity of Theorem 1 are given in subsection C and subsection D, respectively.

A. Tools: Quantum Condition Number Comparison Algorithm
Realizing that multicollinearity appears with large κ [52,53], and hence small eigenvalues, the following preselection algorithm is developed: repeat a simplified phase estimation sub-algorithm until an eigenvalue small enough is detected.If such an eigenvalue is found, the corresponding portfolios will be recorded as an alternative one.It worth noticing that some cointegrated pairs may be missed in our algorithm, but it does not matter since our task is to search for some collinear portfolios instead of the impossible mission to find all of the cointegrated pairs.We denote this procedure Quantum Condition Number Comparator QCN C(κ, ϕ) and get the following result: Lemma 2 Supposing that A is an N × N Hermitian matrix with A = 1 with unknown condition number κ and the probability density function of eigenvalues is p(λ).Then there is a quantum algorithm using ) calls of A to determine whether the condition number is larger than κ 0 .In the case of a uniform probability distribution, A's calls are O(κ 2 0 log (1/ǫ)) so that whenever κ ≥ 2κ 0 , the target qubit will be 1.
It should be noticed that this repeating time, especially when κ is large, is determined by the threshold κ 0 , while traditional algorithms depend on the unknown κ.This is an algorithm finding whether the condition number of a linear system is large than the given threshold without solving the equations.P roof of Lemma 2. Without loss of generality, suppose that A is a matrix with Frobenius norm , and unknown rank r.A direct calculation shows that: = d/r (4) ≥1. ( Here in (2) is the induced L 2 norm and equals to |λ max (A)| (see [54]'s example 5.6.6), and it follows the inequality (3) (see [55]).
For any given eigenvector |λ , with a variant of phase estimation, it is easy for us to determine whether its corresponding eigenvalue λ is larger than 1/κ 0 or not with complexity O(κ 0 log(1/ǫ)) [46].By the definition of the condition number of normal matrices, for any known eigenvalue λ: Hence a lower bound of the condition number is also given.Whenever a sufficiently small eigenvalue λ 0 is given, the matrix can be regarded with condition number greater than κ 0 and high multicollinearity as a consequence.
Obviously, there is a certain probability of success when the testing eigenvalue is larger than κ 0 .Let the condition number be κ and the probability density function of eigenvalues be p; the success probability is Under the assumption that the eigenvalues follow a uniform probability distribution, the success probability turns to be Here κ is assumed large and κ − 1 ≃ κ.Moreover, whenever a matrix with κ ≥ 2κ 0 is given, we have: This procedure shall be repeated 2κ 0 times to boost the success probability.Hence the total number of calls for A is O(κ 2 0 log (1/ǫ)).Since complexity to simulate It should be mentioned that the assumption of uniform distribution is reasonable.Although different distributions of eigenvalues may appear in specified realistic problems, some normal conditions can be imposed to guarantee that the algorithm will still work with slight modification.

B. Algorithm
To see how to derive an algorithm more efficient on κ, one should notice that matrices with small condition numbers can be found quite early and need not be calculated anymore.Some clock registers are used to obsolete those matrices with small condition numbers by starting from small threshold κ 0 .The larger the threshold κ 0 is, the fewer the matrices need to be tested.Repeating this procedure several times can make the acceleration.
Suppose that M = ⌈log κ 0 ⌉.The clock registers C 1,...,M are used to control and store the result of subprocedure A j defined later.Another 1-qubit register F is used as a flag register to donate if the algorithm is stopped.For all j ∈ {1, ..., M }, let φ j = 1/κ j = 2 −j , and let ǫ be the desired precision.Since it is our target to verify whether matrix A contains components corresponding to small eigenvalues, the algorithm is defined as A = A M A M−1 ...A 1 , where A j is defined as follows: Algorithm A j Conditional on first j−1 qubits of H C being |1 , apply QCNC(κ j , ǫ) using C j as the output qubit and additional fresh qubits from P as ancilla (denoted by P j ).If C j is left |0 in the first term, the qubit on stop flag register F will be flipped.

C. Correctness
We shall now prove the correctness of this algorithm.
P roof of Theorem 1 (correctness part).Given a matrix A, the condition number is either in some interval [κ j , 2κ j ] or greater than κ 0 .i) Suppose that a matrix with condition number κ j ≤ κ ≤ 2κ j is given: State after A 1 to A j−1 Since κ ≤ φ 1,...,j−1 , the clock registers C 1 , ..., C j−1 are at position |1 while the stop flag register F stays |0 with high probability.After j-1 steps the state is left as where γ i 1 is the ancillary state produced by the i th call to QCNC.
State after A j Because φ j ≤ κ ≤ 2φ j , QCNC will split the j th control register to |1 Cj with high probability: where U j = 1 j 0 m−j .
Since C j is left |0 in the first term, the qubit on register F is flipped: State after A j+1 This will affect the two parts of the state in different way: In the case that the j th control qubit is splitted, the step QCNC(κ j+1 , ǫ) is implemented.Notice that κ ≤ φ j+1 , the state turns to be: As for the case that j th control qubit is |0 , nothing will be done and the state is:

State after A
Given a matrix A with condition number κ j ≤ κ < κ j+1 , the final state at the end of algorithm A is: ii)As for matrix A with condition number κ ≥ κ 0 , we have: It should be noticed that whenever the flag register is splitted to |1 F , the stops at some step and reject the hypothesis of the condition number larger than κ 0 .a measurement can be implemented on |1 F to decide whether the matrix contains a component with eigenvalues less than some given 1/κ 0 .Besides, a control counter circuit can be employed on the M clock registers to probe the range for κ.

D. Complexity Analysis
In this subsection, the algorithm's complexity analysis is given to finish the proof of Theorem 2. It should be mentioned that the algorithm complexity depends on the specific distribution of condition numbers and eigenvalues of the problem.In this work, a theoretical framework is developed for analysis.Moreover, as an example, the result assuming that log κ follows a uniform probability distribution is calculated.This assumption is common and reasonable since there are relatively fewer matrices with large condition number.P roof of Theorem 1 (complexity part).Suppose that there are n matrices and the condition number threshold for comparison is κ 0 .Let M = ⌈log κ 0 ⌉.Let κ j = 2 j and P j be the probability that the matrix's condition number satisfies κ j−1 ≤ κ ≤ κ j .Then the cummulative number of queries T j for this kind of matrix is: Hence cnosidering the probability, the arithmatic average number of queries is: Supposing that log κ follows a uniform probability contribution, the probability is and the average time is: Hence the complexity is O( √ dκ 2 0 log (1/ǫ) 2 ) as claimed in Theorem 1.
The complexity of this algorithm also depends on the fixed threshold κ 0 instead of an unknown κ.Hence besides the acceleration compared to classical algorithms, the stability and robustness are also improved to satisfy financial problems.

V. QUANTUM COINTEGRATION TEST
To finish the last peice of quantum statistical arbitrage, it needs to be verified whether the preselected matrices contain a cointegrated pair.The global structure and details of QCT are described in the first subsection, and the analysis of complexity is given in the second subsection.These two parts yield the following result: Theorem 3 Suppose that d and N are the number of kinds of stocks and the time length of stock prices, ǫ is the precision desired, and κ is the condition number.Then the cointegration test with L lag-length augmented dickey fuller test can be implemented with complexity O( d 2.5 κ 3 δ 2 poly(log ), where δ = min{1/d, ǫ}, δ ′ = min{1/(L + 2), √ L + 2ǫ 2 }.

A. Algorithm
First of all, the following procedure is used to generate the residual sequence of linear regression.Since the residuals sequence is needed instead of regression coefficients or predicted values [42,43], known quantum linear algorithms should be employed with some further modification.The work of [43]'s Theorem2 is used to derive an approximation β of the regression coefficients β.Lemma 4 (QLR, Theorem 2 in [43]) Let X = (x i,j ) be an N * d balanced matrix such that its singular values are in range [1/κ, 1].Let y = (y 1 , y 2 , ..., y N ) T be a balanced unit vector.Suppose (X, y) is well behaved.Given ǫ > 0 and access to the procedures P x and P y described above.Then the problem to output a vector β = (β 1 , β 2 , ..., β d ) T such that β − β ≤ ǫ and β = X † y can be solved by a gate-efficient quantum algorithm that makes O( d 2.5 κ 3 δ 2 poly[log 2 ( dκ δ )]) uses of P x and P y , where δ = min{1/d, ǫ}.This Quantum Linear Regression procedure is denoted as QLR(d, δ, κ).Then the predicted value vector ŷ is calculated by the matrix multiplication ŷ = Xβ, and the residuals sequence is derived by a vector subtraction between the predicted values ŷ and real values y: This should be a hybrid algorithm since classical algorithms can calculate matrix multiplications and subtractions with fewer restrictions and more efficiently.Next, another regression QLR(L + 1, δ ′ , κ ′ ) on time variable and lagged residuals will be employed to derive the statistical index.The lagged residuals ∆u t is defined as the first-order difference and can be calculated efficiently by a vector subtraction: Then QLR(L + 1, δ ′ , κ ′ ) procedure shows: where L is the lag-length used in the ADF test, and β is the coefficient of the time variable t.The test statistic DF T = γ SE(γ) , where SE means standard error, can be computed by classcial computer more efficiently.
Finally, the result will be sent to be compared with a critical value table [38].And the total algorithm is summarised as Algorithm 3.

Input:
κ 0 : the threshold for preselection T : the length of time interval J: the total number of stocks p (j) t : the j th stock's price at time t.Output: (f, β)flag and cointegrated coefficients.
Data Loading: By (41), the classical benchmark is O(N 2 d) = 10 18 , and the average complexity of our algorithm is mainly determined by the first preselection subroutine wtih complexity O( √ dκ 2 0 log (1/ǫ) 2 )) = 10 8 (see details below).The primary reason for this acceleration is that the preselection procedure can search the multicollinearity without large matrix factorizations and regressions.Secondly, the problem size determined by sample number N is supposed to be very large for our problem of high-frequency trading: On the one hand, for high-frequency trading, there is a short time interval and a large number N 0 of trading date quotes of every single trading day.On the other hand, it does make sense in finance to consider a long time interval l since it is a statistical arbitrage model instead of some models for prediction such as momentum trading.Finally, for the specific case of statistical arbitrage trading strategy, it is common and unavoidable to handle matrices with large condition number κ, resulting in high cost of computing resources and time complexity.Utilizing the ability to detect κ by QCNCA, our algorithm is time variable one and adaptive to κ.Since most portfolios are with small κ as discussed above, giving a bound κ 0 = 1000, our algorithm's complexity is about O( √ dκ 2 0 log (1/ǫ) 2 )) = 10 8 .
The number of qubits needed can be estimated as follows: According to the data size discussed above, the qubits needed to prepare for the initial state is about log (1.2 * 10 7 ) + log 8000 ≈ 35.The QCN C(κ, φ) circuit consists of simplified phase estimation subcircuits, and each subcircuit with 0.1 precision needs more than 4 qubits.Moreover, the V T P A circuit consists of QCN C(κ j , φ) circuits for different κ j , and hence more than 50 qubits are needed, which are hard for us to simulate.

VII. CONCLUSION
In this article, we introduce quantum algorithms for quantitative trading in the case of high-frequency statistical arbitrage and show the quantum advantage.Besides wxploring new financial applications, two heuristic algorithms are also developed as instruments: One is for the estimation of the condition number of a given matrix, which has not been considered and proposed before as far as we know.This algorithm can be applied to solve other problems where condition number is a primary influencing factor of the algorithm's complexity, such as quantum computational fluid dynamics and differential equation solution [60][61][62][63][64].The other is the implemention of statistical cointegration test, which has many applications in time series, finance analysis.Some modifications and exploration will be considered later to suit these exciting problems.
During the analysis of QCNCA and VTPA's complexity, we provide a theoretical framework and show the quantum advantage under the assumption of uniform distribution.Since the real problems are complicated, many other statistical models and different distributions will be taken into consideration.By some modification in Eqs.(8)(9)(10)(11), this method might still work with different results of complexity, and this is our further research direction.Moreover, the work of circuit simplification and simulation will be done in the future.