Option Pricing using Quantum Computers

We present a methodology to price options and portfolios of options on a gate-based quantum computer using amplitude estimation, an algorithm which provides a quadratic speedup compared to classical Monte Carlo methods. The options that we cover include vanilla options, multi-asset options and path-dependent options such as barrier options. We put an emphasis on the implementation of the quantum circuits required to build the input states and operators needed by amplitude estimation to price the different option types. Additionally, we show simulation results to highlight how the circuits that we implement price the different option contracts. Finally, we examine the performance of option pricing circuits on quantum hardware using the IBM Q Tokyo quantum device. We employ a simple, yet effective, error mitigation scheme that allows us to significantly reduce the errors arising from noisy two-qubit gates.


INTRODUCTION
Options are financial derivative contracts that give the buyer the right, but not the obligation, to buy (call option) or sell (put option) an underlying asset at an agreed-upon price (strike) and timeframe (exercise window). In their simplest form, the strike price is a fixed value and the timeframe is a single point in time, but exotic variants may be defined on more than one underlying asset, the strike price can be a function of several market parameters and could allow for multiple exercise dates. As well as providing investors with a vehicle to profit by taking a view on the market or exploit arbitrage opportunities, options are core to various hedging strategies and as such, understanding their properties is a fundamental objective of financial engineering. For an overview of option types, features and uses, we refer the reader to Ref. [1].
Due to the stochastic nature of the parameters options are defined on, calculating their fair value can be an arduous task and while analytical models exist for the simplest types of options [2], the simplifying assumptions on the market dynamics required for the models to provide closed-form solutions often limit their applicability [3]. Hence, more often than not, numerical methods have to be employed for option pricing, with Monte Carlo being one of the most popular due to its flexibility and ability to generically handle stochastic parameters [4,5]. However, despite their attractive features in option pricing, classical Monte Carlo methods generally require extensive computational resources to provide accurate option price estimates, particularly for complex options. Because of the widespread use of options in the finance industry, accelerating their convergence can have a significant impact in the operations of a financial institution.
By leveraging the laws of quantum mechanics a quantum computer [6] may provide novel ways to solve computationally intensive problems such as quantum chemistry [7][8][9][10], solving linear systems of equations [11], and machine learning [12][13][14]. Quantitative finance, a field with many computationally hard problems, may benefit from quantum computing. Recently developed applications of gate-based quantum computing for use in finance [15] include portfolio optimization [16], the calculation of risk measures [17] and pricing derivatives [18][19][20]. Several of these applications are based on the Amplitude Estimation algorithm [21] which can estimate a parameter with a convergence rate of 1/M , where M is the number of quantum samples used. This represents a theoretical quadratic speed-up compared to classical Monte Carlo methods.
In this paper we extend the pricing methodology presented in [17,18] and place a strong emphasis on the implementation of the algorithms in a gate-based quantum computer. We first classify options according to their features and show how to take the different features into account in a quantum computing setting. In Section III, we review the quantum algorithms needed to price options and discuss how to represent relevant probability distributions in a quantum computer. In Section IV, we show a framework to price vanilla options and portfolios of vanilla options, options with path-dependent payoffs and options on several underlying assets. In Section V we show results from evaluating our option circuits on quantum hardware, and describe the error mitigation scheme we employ to increase the accuracy of the estimated option prices. In particular, we employ the maximum likelihood estimation method introduced in [22] to perform amplitude estimation without phase estimation in option pricing using three qubits of a real quantum device.

II. REVIEW OF OPTION TYPES AND THEIR CHALLENGES
We classify options according to two categories: pathindependent vs path-dependent and options on single assets or on multiple assets. Path-independent options have a payoff function that depends on an underlying asset at a single point in time. Therefore, the price of the asset arXiv:1905.02666v2 [quant-ph] 4 Jul 2019 up to the exercise date of the option is irrelevant for the option price. By contrast, the payoff of path-dependent options depends on the evolution of the price of the asset and its history up to the exercise date. Table I exemplifies this classification. Options that are path-independent and rely on a single asset are the easiest to price. This is done using Amplitude Estimation once a proper representation of the distribution of the underlying asset can be loaded to the quantum computer. Path-independent options on multiple assets are only slightly harder to price since more than one asset is now involved and the probability distribution loaded into the quantum computer must accout for correlations between the assets. Path-dependent options are harder to price than pathindependent options since they require a representation of the possible paths the underlying assets can take in the quantum computer.

III. IMPLEMENTATION ON A GATE BASED QUANTUM COMPUTER
Here we review some of the building blocks needed to price options on a gate-based quantum computer.

A. Distribution loading
The analytical formulas used to price options in the Black-Scholes-Merton (BSM) model [2,23] assume that the underlying stock prices at maturity follow a lognormal distribution with constant volatility. Such distributions can be efficiently loaded in a gate-based quantum computer [18,24]. However, to properly model the market prices of options, the volatility of the geometric brownian process describing the evolution of the assets must be changed for options with different strike prices [25]. This discrepancy between the BSM model and market prices is because stocks do not follow a geometric Browninan motion process with constant volatility. It is thus important to be able to efficiently represent arbitrary distributions of financial data in a quantum computer.
The loading of arbitrary states into quantum systems requires exponentially many gates [26], making it inefficient to model arbitrary distributions as quantum gates. Since the distributions of interest are often of a special form, the limitation may be overcome by using quantum Generative Adverserial Networks (qGAN). These networks allow us to load a distribution using a polynomial number of gates [19]. A qGAN can learn the random distribution X underlying the observed data samples { x 0 , . . . , x k−1 } and load it directly into a quantum state. This generative model employs the interplay of a classical discriminator, a neural network [27], and a quantum generator (a parametrized quantum circuit). More specifically, the qGAN training consists of alternating optimization steps of the discriminator's parameters φ and the generator's parameters θ. After the training, the output of the generator is a quantum state that represents the target distribution. The n-qubit state |i n = |i n−1 ...i 0 encodes the integer i = 2 n−1 i n−1 + ... + 2i 1 + i 0 ∈ {0, ..., 2 n − 1} with i k ∈ {0, 1} and k = 0, ..., n − 1. The probabilities p i (θ) approximate the random distribution underlying the training data. We note that the outcomes of a random variable X can be mapped to the integer set {0, ..., 2 n − 1} using an affine mapping. Furthermore, the approach can be easily extend to multivariate data, where we use a separate register of qubits for each dimension [19].

B. Amplitude Estimation
The advantage of pricing options on a quantum computer comes from the amplitude estimation (AE) algorithm [21] which provides a quadratic speed-up over classical Monte-Carlo simulations [28,29]. Suppose a unitary operator A acting on a register of (n+1) qubits such that for some normalized states |ψ 0 n and |ψ 1 n , where a ∈ [0, 1] is unknown. AE allows the efficient estimation of a, i.e., the probability of measuring |1 in the last qubit. This estimation is obtained with an operator Q, based on A, and Quantum Phase Estimation [30] to approximate FIG. 2. Quantum circuit that creates the state in Eq. (4). Here, the independent variable i = 4i2 + 2i1 + i0 ∈ {0, ..., 7} is encoded by three qubits in the state |i 3 = |i2i1i0 with i k ∈ {0, 1}. Therefore, the linear function f (i) = f1i + f0 is given by 4f1i2 + 2f1i1 + f1i0 + f0. After applying this circuit the quantum state is |i 3 The circuit on the right shows an abbreviated notation.
with probability of at least 8/π 2 . This represents a quadratic speedup compared to the O M −1/2 convergence rate of classical Monte Carlo methods [31].

C. Linearly controlled Y-rotations
We obtain the expectation value of a linear function f of a random variable X with AE by creating the operator A such that a = E[f (X)], see Eq. (2). Once A is implemented we can prepare the state in Eq. (2) and the Q operator. In this section, we show how to create a close relative of the operator in Eq. (2) and then, in Section III D, we show how to use AE.
Since the payoff function for option portfolios is piecewise linear we only need to consider linear functions f : We can efficiently create an operator that performs using controlled Y-rotations [17]. To implement the linear term of f (i) each qubit j (where j ∈ {0, . . . n − 1}) in the |i n register acts as a control for a Y-rotation with angle 2 j f 1 of the ancilla qubit. The constant term f 0 is implemented by a rotation of the ancilla qubit without any controls, see Fig. 2. The controlled Y-rotations can be implemented with CNOT and single-qubit gates [32].

D. Expectation value of functions using AE
We now describe how to obtain E[f (X)] for a linear function f of a random variable X which is mapped to integer values i ∈ {0, ..., 2 n − 1} that occur with probability p i . To do this we create the operator that maps i √ p i |i n |0 to using the procedure outlined in Sec. III C. The parameter c ∈ [0, 1] is a scaling parameter. The functionsf (i) and f (i) are related bỹ Here f min = min i f (i) and f max = max i f (i). The relation in Eq. (5) is chosen so thatf (i) ∈ [−1, 1]. Thus, sin 2 [cf (i) + π/4] is an anti-symmetric function around f (i) = 0. With these definitions, the probability to find the ancilla qubit in state |1 , namely is well approximated by To obtain this result we made use of the approximation which is valid for small values of cf (i). With this first order approximation the convergence rate of AE is O(M −2/3 ) when c is properly chosen which is already faster than classical Monte Carlo methods [17]. We can recover the O(M −1 ) convergence rate of AE by using higher orders implemented with quantum arithmetic. The resulting circuits, however, have more gates. This trade-off, discussed in Ref. [17], also gives a formula that specifies which value of c to use to minimize the estimation error made when using AE. From Eq. (6) we can recover E[f (X)] since AE allows us to efficiently retrieve P 1 and because we know the values of f min , f max and c.

IV. OPTION PRICING ON A QUANTUM COMPUTER
In this section we show how to price the different options shown in Tab. I. We put an emphasis on the implementation of the quantum circuits that prepare the states needed by AE. We use the different building blocks reviewed in Sec. III.

A. Path-independent options
The price of path-independent vanilla options (e.g. European call and put options) depend only on the distribution of the underlying asset price S T at the option maturity T and the payoff function f (S T ) of the option. To encode the distribution of S T in a quantum state we truncate it to the range [S T,min , S T,max ] and discretize this interval to {0, ..., 2 n − 1} using n qubits. In the quantum computer the distribution loading operator P X creates a state with i ∈ {0, ..., 2 n − 1} to represent S T . This state, exemplified in Fig. 3, may be created using the methods discussed in Sec. III A. We start by showing how to price vanilla call or put options and then generalize our method to capture the payoff structure of portfolios containing more than one vanilla option. FIG. 3. Example price distribution at maturity loaded in a three-qubit register. In this example we followed the Black-Scholes-Merton model which implies a lognormal distribution of the asset price at maturity T with probability density func- . σ is the volatility of the asset and µ = r − 0.5σ 2 T + ln(S0), with r the riskfree market rate and S0 the asset's spot at t = 0. In this figure we used S0 = 2, σ = 10%, r = 4% and T = 300/365.

Vanilla options
To price vanilla options with strike K, we implement a comparison between the values in state (8) with K. A quantum comparator circuit sets an ancilla qubit |c , initially in state |0 , to the state |1 if i ≥ K and |0 otherwise. The state |ψ n in the quantum computer therefore undergoes the transformation This operation can be implemented by a quantum comparator [33] based on CNOT and Toffoli gates. Since we know the value of the strike, we can implement a circuit tailored to the specific strike price. We use n ancilla qubits |a 1 , ..., a n and compute the two's complement of the strike price K in binary using n bits, storing the digits in a (classical) array t[n]. For each qubit |i k in the |i n register, with k ∈ {0, ..., n − 1}, we compute the possible carry bit of the bitwise addition of |i k and t[k] into |a k . If t[k] = 0, there is a carry qubit at position k only if there is a carry at position k − 1 and |i k = 1. If t[k] = 1, there is a carry qubit at position k if there is a carry at position k − 1 or |i k = 1. After going through all n qubits from least to most significant, |i n will be greater or equal to the strike price, only if there is a carry at the last (most significant) qubit. This procedure along with the necessary gate operations is illustrated in Fig. 4. An implementation for K = 1.9 and a three-qubit register is shown in Fig. 6.
To prepare the operator for use with AE we add to |φ 1 a second ancilla qubit initially in the state cos(g 0 ) |0 + sin(g 0 ) |1 . Here, g 0 is an angle with a value that we will carefully select. Next, we perform a rotation of the new ancilla qubit controlled by the comparator qubit |c and the qubits in |ψ n . The state This operation, implemented by the quantum circuit in Fig. 7, applies a rotation with an angle g(i) only if i ≥ K. The probability to find the second ancilla in state |1 , efficiently measurable using AE, is Now, we must carefully chose the angle g 0 and the function g(i) to recover the expected payoff E[f (X)] of the option from P 1 using the approximation in Eq. (6). The payoff function of vanilla options is piece-wise linear We now focus on a European call option with payoff Sec. III D, we must set Circuit that compares the value represented by an n-qubit register |i n , to a fixed value K. We use n ancilla qubits |a1, ..., an , a classical array t[n] holding the precomputed binary value of K's two's complement and a qubit |c which will hold the result of the comparison with |c = 1 if |i ≥ K. For each qubit |i k , with k ∈ {1, ..., n}, we use a Toffoli gate to compute the carry at position k if t[k] = 1 and a logical OR, see Circuit that computes the logical OR between qubits |a and |b into qubit |c . The circuit on the right shows the abbreviated notation used in Fig. 4. 6. Quantum circuit that sets a comparator qubit |c to |1 if the value represented by |i 3 is larger than a strike K = 1.9, for the spot distribution in Fig. 3. The unitary PX represents the set of gates that load the probability distribution in Eq. (8). An ancilla qubit |a is needed to perform the comparison. It is uncomputed at the end of the circuit.
where i max = 2 n − 1. This choice of g(i) forces us to Circuit that creates the state in Eq. (9). We apply this circuit directly after the comparator circuit shown in Fig. 6. The multi-controlled y-rotation is the gate shown in Fig. 2 controlled by the ancilla qubit |c that contains the result of the comparison between i and K. chose To see why, we substitute Eqs. (12) and (13) in Eq. (10) and use the approximation in Eq. (7). Therefore, This shows us that we needed g 0 = π/4 − c to used the identity up to a scaling factor and a constant. From this last equality we recover the expected payoff of the option given the probability distribution of the underlying asset. We should note that the fair value of the option requires appropriately discounting the expected payoff of the option to today, but as the discounting can be performed after the expectation value has been calculated we omit it from our discussion for simplicity. We demonstrate the performance of our approach by running amplitude estimation using Qiskit [34] on the overall circuit produced by the elements described in this section, and verifying the convergence to the analytically computed value or classical Monte Carlo estimate. An illustration of the convergence of a European call option with increasing evaluation qubits is shown in Fig. 8.
A straightforward extension of the analysis above yields a pricing model for a European put option, whose

Portfolios of options
Various popular trading and hedging strategies rely on entering multiple option contracts at the same time in-stead of individual call or put options and as such, these strategies allow an investor to effectively construct a payoff that is more complex than that of vanilla options. For example, an investor that wants to profit from a volatile asset without picking a direction of where the volatility may drive the asset's price, may choose to enter a straddle option strategy, by buying both a call and a put option on the asset with the same expiration date and strike. If the underlying asset moves sharply up to expiration date, the investor can make a profit regardless of whether it moves higher or lower in value. Alternatively, the investor may opt for a butterfly option strategy by entering four appropriately structured option contracts with different strikes simultaneously. Because these option strategies give rise to piecewise linear payoff functions, the methodology described in the previous section can be extended to calculate the fair values of these option portfolios.
In order to capture the structure of such option strategies, we can think of the individual options as defining one or more effective strike prices K j and add a linear function f j (S) = a j S + b j between each of these strikes. For example, to price an option strategy with the payoff function which corresponds to a call spread (the option holder has purchased a call with strike K 1 and sold a call with strike K 2 ), we use the functions f 0 , f 1 , and f 2 such that To match Eq. (15) with Eq. (16) we set f 0 (S) = 0, f 1 (S) = S − K 1 and f 2 (S) = −S + K 2 . In general, to price a portfolio of options with m effective-strike prices K 1 , ..., K m and m+1 functions f 0 (S), ..., f m (S) we need an ancilla qubit per strike to indicate if the underlying has reached the strike. This allows us to generalize the discussion from Sec. IV A 1. We apply a multi-controlled Y-rotation with angle g j (i) if i ≥ K j for each strike K j with j ∈ {1, ..., m}. The rotation g 0 (i) is always applied, see the circuit in Fig. 9. The functions g j (i) are determined using the same procedure as in Sec. IV A 1.

B. Multi-asset and path-dependent options
In this section we show how to price options with pathdependent payoffs as well as options on more than one underlying asset. In these cases, the payoff function depends on a multivariate distribution of random variables {S j } with j ∈ {1, ..., d}. The S j 's may represent one or several assets at discrete moments in time or a basket of assets at the option maturity. In both cases, the probability distribution of the random variables S j are truncated to the interval [S j,min , S j,max ] and discretized using Quantum circuit that implements the multicontrolled Y-rotations for a portfolio of options with m strike prices.
2 nj points so that they can be represented by d quantum registers where register j has n j qubits. Thus, the multivariate distribution is represented by the probabilities p i1,...,i d that the underlying has taken the values i 1 , ..., i d with i j ∈ {0, ..., 2 nj − 1}. The quantum state that represents this probability distribution, a generalization of Eq. (8), is with n = j n j . Various types of options, such as Asian options or basket options, require us to compute the sum of the random variables S j . The addition of the values in two quantum registers |a, b → |a, a + b may be calculated in quantum computers with adder circuits based on CNOT and Toffoli gates [35][36][37]. To this end we add an extra qubit register with n qubits to serve as an accumulator. By recursively applying adder circuits we perform the transformation |ψ n |0 n → |φ n+n where |φ n+n is given by Here circuit optimization may allow us to perform this computation in-place to minimize the number of qubit registers needed. Now, we use the methods discussed in the previous section to encode the option payoffs into the quantum circuit.

Basket Options
A European style basket option is an extension of the single asset European option discussed in Sec. IV A, only now the payoff depends on a weighted sum of d underlying assets. A call option on a basket has the payoff profile where S basket = w · S, for basket weights w = [w 1 , w 2 , . . . , w d ], w i ∈ [0, 1], underlying asset prices at option maturity S = [S 1 , S 2 , . . . S d ] and strike K. In the BSM model, the underlying asset prices are described by a multivariate lognormal distribution with probability density function [38] where ln S = (ln S 1 , ln S 2 . . . , ln S d ) T and µ = (µ 1 , µ 2 . . . µ d ) T , where each µ i is the lognormal distribution parameter for each asset defined in the caption of Fig. 3. Σ is the d × d positive-definite covariance matrix of the d underlyings with σ i the volatility of the ith asset, −1 ≤ ρ ij ≤ 1 the correlation between assets i and j and T the time to maturity. The quantum circuit for pricing a European style basket call option is analogous to the single asset case, with an additional unitary to compute the weighted sum of the uncertainty registers |i 1 n1 . . . |i d n d before applying the comparator and payoff circuits, controlled by the accumulator register |b n = |i 1 + ... + i d n . A schematic of these components is shown in Fig. 10. The implementation details of the circuit that performs the weighted sum operator is discussed in Appendix A.
We use a basket option to compare the estimation accuracy between AE and classical Monte Carlo. For M applications of the Q operator, see Fig. 1, the possible values returned by AE will be of the form sin 2 (yπ/M ) for y ∈ {0, ..., M − 1} and the maximum distance between two consecutive values is This quantity determines how close M operations of Q can get us to the amplitude which encodes the option price. Using sin 2 (π/4+x) = x+1/2+O( (3) and Eq. (14) we can determine that with probability of at least 8/π 2 , our estimated option price using AE will be within of the exact option price, where c, i max and K are the parameters used to encode the option payoff into our quantum circuit, discussed in Sec. IV A 1. To compare this estimation error to Monte Carlo, we use the same number of samples to price an option classically, and determine the approximation error at the same 8/π 2 ≈ 81% confidence interval by repeated simulations. The comparison for a basket option on three underlying assets shows that AE provides a quadratic speed-up, see Fig. 11.

Asian Options
We now examine arithmetic average Asian options which are single-asset, path-dependent options whose payoff depends on the price of the underlying asset at multiple time points before the option's expiration date. Specifically, the payoff of an Asian call option is given by where K is the strike price,S is the arithmetic average of the asset's value over a pre-defined number of points d between 0 and the option maturity T The probability distribution of asset prices at time t will again be lognormal with probability density function with µ t = r − 0.5σ 2 ∆t + ln(S t−1 ) and ∆t = T /d. We can then use the multivariate distribution in Eq. (20), with S now a d-dimensional vector of asset prices at time points [t 1 . . . t d ], instead of distinct underlying prices at maturity T . As we are not considering multiple underlying assets that could be correlated, the covariance matrix is diagonal Σ = ∆t[diag(σ 2 , ..., σ 2 )]. An illustration of the probability density function used for an asset defined on two time steps is shown in Fig. 12.
We now prepare the state |ψ n , see Eq. (17), where each register represents the asset price at each time step up to maturity. Using the weighted sum operator of Appendix A with equal weights 1/d, we then calculate the average value of the asset until maturity T , see Eq. (25), into a register |S Finally, we use the same comparator and rotation circuits that we employed for the basket option illustrated in Fig. 10 to load the payoff of an arithmetic average Asian option into the payoff qubit |p .

Barrier Options
Barrier options are another class of popular option types whose payoff is similar to vanilla European Op- • Knock-In: The option has no value unless the underlying asset crosses a certain price level before maturity.
If the required barrier event for the option to have value at maturity occurs, the payoff then depends only on the value of the underlying asset at maturity and not on the path of the asset until then. If we consider a Knock-In barrier option and label the barrier level B, we can write the option's payoff as where T is the time to maturity, S t the asset price at time t with 0 < t ≤ T and K the option strike.
To construct a quantum circuit to price a Knock-In barrier option, we use the same method as for the Asian option where T is divided into d equidistant time intervals with ∆t = T /d, and use registers |i 1 n1 |i 2 n2 . . . |i d n d to represent the discretized range of asset prices at time t ∈ {∆t, 2∆t, . . . , d·∆t = T }. The probability distribution of Eq. (26) is used again to create the state |ψ n in Eq. (17).
To capture the path dependence introduced by the barrier, we use an additional d-qubit register |b d to monitor if the barrier is crossed. Each qubit |b t in |b d is set to |1 if |i t nt ≥ B. An ancilla qubit |b | is set to |1 if the barrier has been crossed in at least one time step. This is done by computing the logical OR, see Fig. 5, of every qubit in |b d and storing the result in the ancilla This is computed with X (NOT) and Toffoli gates and d − 2 ancilla qubits. The ancilla qubit |b | is then used as a control for the payoff rotation into the payoff qubit, effectively knocking the option in. For Knock-Out barrier options, we can follow the same steps and apply an X gate to the ancilla barrier qubit before using it as control, in this manner knocking the option out if the barrier level has been crossed. A circuit displaying all the components required to price a Knock-In barrier option is shown in Fig. 13. Results from amplitude estimation on a barrier option circuit using a quantum simulator are shown in Fig. 14.
Even though we have focused our attention on barrier options where the barrier event is the underlying asset crossing a barrier from below, we can use the same method to price barrier options where barrier events are defined as the asset crossing the value from above. This only requires changing the comparator circuits to compute S t ≤ B in the barrier register |b d .

V. QUANTUM HARDWARE RESULTS
In this section we show results for a European call option evaluated on quantum hardware. We use three qubits, two of which represent the uncertainty and one encodes the payoff.
To examine the behavior of the circuit for different input probability distributions, we run eight experiments that differ by the initial spot price S 0 and all other parameters are kept constant. The spot price is varied from 1.8 to 2.5 in increments of 0.1. This way we can use the same circuit for all experiments and only vary the Y-rotation angles used to map the initial probability distribution onto the qubit register. This choice of input parameters allows us to evaluate our circuits for expected option prices in the range [0.0754, 0.7338].
Each experiment is evaluated on the IBM Q Tokyo 20qubit device with 8192 shots. We repeat each 8192-shot experiment 20 times and average over the 20 measured probabilities in order to limit the effect of any one-off issues with the device. The standard deviation of the measured probabilities across the 20 runs was always < 2%. The connectivity of IBM Q Tokyo allows to choose three fully connected qubits for the experiments, and thus, no swaps are required to implement any two-qubit gate in our circuits [34]. For all circuits described in the following sections, we used qubits 6, 10 and 11.

A. Algorithm and Operators
We now show how to construct the operator A that is required for AE. The log-normal distribution on two qubits can be loaded using a single CNOT gate and four single qubit rotations [39]. To encode the payoff function we also exploit the small number of qubits and apply a uniformly controlled Y-rotation instead of the generic construction using comparators introduced in Sec. IV. A uniformly controlled Y-rotation, i.e.
implements a different rotation angle θ i , i = 0, ..., 2 n − 1 for each state of the n-control qubits. For n = 2, this operation can be efficiently implemented using four CNOT gates and four single qubit Y-rotations [40,41]. The resulting circuit implementing A is shown in Fig. 15. Note that in our setup different initial distributions only lead to different angles of the first four Y-rotations and do not affect the rest of the circuit. Although we use a uniformly controlled rotation, the rotation angles are constructed in the same way described in Sec. III D. We use an approximation scaling of c = 0.25 and the resulting angles are [θ 0 , . . . , θ 3 ] = [1.1781, 1.1781, 1.5708, 1.9635], which shows the piecewise linear structure of the payoff function.
The total resulting circuit requires five CNOT gates and eight single-qubit Y-rotations, see Fig. 15. Since we use uniformly controlled rotations, we do not need any ancilla qubit. Note that if we want to evaluate the circuit for A alone, we can replace the last CNOT gate in Fig. 15 by classical post-processing of the measurement result: if q 1 is measured as |1 , we flip q 2 and otherwise we do nothing. This further reduces the overall CNOT gate count to four.
A quadratic speed-up can also be realized by performing AE without quantum phase estimation [22]. This is done by measuring Q k A |0 for k = 2 0 , ..., 2 m−1 for a given m and applying a maximum likelihood estimation. If we define M = 2 m − 1, i.e. the total number of Q-applications, and we consider N shots for each experi-  15. The A operator of the considered European call option: first, the 2-qubit approximation of a log-normal distribution is loaded, and second, the piecewise linear payoff function is applied to last qubit controlled by the first two. This operator can be used within amplitude estimation to evaluate the expected payoff of the corresponding option. ment, it has been shown that the resulting estimation error scales as O(1/(M √ N )), i.e., the algorithms achieves the quadratic speed-up in terms of M . This leads to shorter circuits than the original implementation of AE, see Appendix B for more details. In the remainder of this section, we focus on QA |0 3 , i.e., the outlined algorithm for m = 1, to demonstrate option pricing on real quantum hardware.
After optimzing the gate count, the resulting circuit for QA |0 3 consists of 18 CNOT gates and 33 single-qubit gates. The detailed circuit diagram and applied circuit optimization steps are provided in Appendix C.

B. Error mitigation and results
We run the circuits for A |0 3 and QA |0 3 on noisy quantum hardware. The results are affected by readout errors and errors that occur during the execution of the circuits.
To mitigate readout errors we run a calibration sequence in which we individually prepare and measure all eight basis states [34,42]. The result is a 8 × 8 readoutmatrix R that holds the probability of measuring each basis state as function of the basis state in which the system was prepared. We use R to correct all subsequent measurements. The error we measure on P 1 for A |0 3 was reduced from ∼ 6% to ∼ 4% using readout error mitigation.
Errors occuring during the quantum circuit can be mitigated using Richardson extrapolation [43]. First, the quantum circuit is run using a rescaled Hamiltonian to amplify the effect of the noise. Second, a Richardson extrapolation is used to extract the result of the quantum circuit at the zero noise limit. In hardware, error mitigation is done by stretching the duration of the gates. For each stretch factor the qubit gates need to be recalibrated [8]. Here, we use a simplified error mitigation protocol that circumvents the need to recalibrate the gates but still allows us to increase the accuracy of the quantum hardware [44]. Since the single-qubit and CNOT gates have an average randomized benchmarking fidelity of 99.7% and 97.8%, respectively, we restrict our error mitigation to the CNOT gates. Furthermore, because the optimized circuit for A |0 3 contains only 4 CNOT gates, we employ the error mitigation protocol only when evaluating QA |0 3 which consists of 18 CNOT gates.
We run the circuit for QA |0 3 three times. In each run we replace the CNOT gates of the original circuit by one, three and five CNOT gates for a total of 18, 54, and 90 CNOT gates, respectively. Since a pair of perfect CNOT gates simplifies to the identity these extra gates allow us to amplify the error of the CNOT gate without having to stretch the gate duration, thus, avoiding the need to recalibrate the gate parameters. As the number of CNOT gates is increased the probability of measuring |1 tends towards 0.5 for all initial spot prices, see Fig. 16(b). After applying the Richardson extrapolation we recover the same behavior as the option price obtained from classical simulations, see Fig. 16(c). Our simple error mitigation scheme dramatically increased the accuracy of the calculated option price: it reduced the error, averaged over the initial spot price, from 62% to 21%.

VI. CONCLUSION
We have presented a methodology and the quantum circuits to price options and option portfolios on a gatebased quantum computer. We showed how to account for some of the more complex features present in exotic options such as path-dependency with barriers and averages. The results that we show are available in Qiskit Finance [34]. Future work may involve calculating the price derivatives [45] with a quantum computer. Pricing options relies on AE. This quantum algorithm allows a quadratic speed-up compared to traditional Monte Carlo simulations but will most likely require a universal faulttolerant quantum computer [46]. However, research to improve the algorithms is ongoing [47][48][49]. Here we have used a new algorithm [22] that retains the AE speedup but that uses less gates to measure the price of an  Fig. 15) (b) Probability of measuring |1 for the QA |0 3 circuit (see Fig. 19). We show the measured probabilities when replacing each CNOT by one, three and five CNOT gates (green, orange, red, respectively), the zero-noise limit calculated using a second-order Richardson extrapolation method (purple), and the probability measured using the simulator (blue). (c) Option price estimated with maximum likelihood estimation from measurements of QA |0 3 and A |0 3 with error mitigation (purple) and without (green). The exact option price for each initial spot price S0 is shown in blue.
option. Furthermore, we employed a simple error mitigation scheme that allowed us to greatly reduce the errors from the noisy quantum hardware. However, larger quantum hardware capable of running deeper quantum circuits with more qubits than the currently available quantum computers is needed to price the typical portfolios seen in the financial industry. Future work could focus on reducing the number of quantum registers in our implementation by performing some of the computation in-place.

VII. ACKNOWLEDGMENTS
The authors want to thank Abhinav Kandala for the very constructive discussions on error mitigation and real quantum hardware experiments.
Opinions and estimates constitute our judgment as of the date of this Material, are for informational purposes only and are subject to change without notice. This Material is not the product of J.P. Morgan's Research Department and therefore, has not been prepared in accordance with legal requirements to promote the independence of research, including but not limited to, the prohibition on the dealing ahead of the dissemination of investment research. This Material is not intended as research, a recommendation, advice, offer or solicitation for the purchase or sale of any financial product or service, or to be used in any way for evaluating the merits of participating in any transaction. It is not a research report and is not intended as such. Past performance is not indicative of future results. Please consult your own advisors regarding legal, tax, accounting or any other aspects including suitability implications for your particular circumstances. J.P. Morgan disclaims any responsibility or liability whatsoever for the quality, accuracy or completeness of the information herein, and for any reliance on, or use of this material in any way. Important disclosures at: www.jpmorgan.com/disclosures IBM, IBM Q, Qiskit are trademarks of International Business Machines Corporation, registered in many jurisdictions worldwide. Other product or service names may be trademarks or service marks of IBM or other companies.
Appendix A: Circuit implementation of weighted sum operator

Weighted sum of single qubits
In this appendix, we demonstrate an implementation of the weighted sum operator on a quantum circuit. The weighted sum operator S computes the arithmetic sum of the values of n qubits |a n = |a 1 . . . a n weighted by n classically defined non-negative integer weights ω = (ω 1 , ω 2 , . . . , ω n ), and stores the result into another mqubit register |s m = |s 1 · · · s m initialized to |0 m . In other words, where The choice of m ensures that the sum register |s m is large enough to hold the largest possible weighted sum, i.e. the sum of all weights. Alternatively, we can write the weights in the form of a binary matrix Ω = (Ω i,j ) ∈ {0, 1} n×n * , where the i-th row in Ω is the binary representation of weight ω i and n * = max d i=1 n i . We use the convention that less significant digits have smaller indices, so |s 1 and Ω i,1 are the least significant digits of the respective binary numbers. Using this binary matrix representation, S is to add the i-th qubit |a i of the state register to the j-th qubit |s j of the sum register if and only if Ω i,j = 1. Depending on the values of the weights, an additional quantum register may be necessary to temporarily store the carries during addition operations. We use |c j to denote the ancilla qubit used to store the carry from adding a digit to |s j . These ancilla qubits are initialized to |0 and will be reset to their initial states at the end of the computation. Based on the above setup, we build quantum circuits for the weighted sum operator using three elementary gates: X (NOT), CNOT, and the Toffoli gate (CCNOT). These three gates suffice to build any Boolean function [35]. Starting from the first column in Ω, for each column j, we find all elements with Ω i,j = 1 and add the corresponding state qubit |a i to |s j . The addition of two qubits involves three operations detailed in Fig. 17: (a) computation of the carry using a Toffoli gate (M), (b) computation of the current digit using a CNOT (D), (c) reset of the carry computation using two X gates and one Toffoli gate ( M). When adding |a i to the j-th qubit of the sum register, the computation starts by applying M and then D to |a i , |s j and |c j , which adds |a i to |s j and stores the carry into |c j . Then, using the same two operations, it adds the carry |c j to the next sum qubit |s j+1 with carry recorded in |c j+1 . The process is iterated until all carries are handled. Finally, it resets the carry qubits by applying M in reverse order of the carry computation. We reset the carry qubits in order to reuse them in later computations if necessary.
In general, we need max(k − 2, 0) carry qubits to com-pute the addition of |a i on |s j , where k ≥ 1 is the smallest integer satisfying k 1|ρ s j,j+k−1 |1 k = 0, where ρ s j,j+k−1 is the density matrix corresponding to |s j · · · s j+k−1 . In the k = 1 case, i.e. |s j = 0, the computation is reduced to "copying" |a i to |s j using the bit addition operator D, and no carries would be produced. For k ≥ 2, Eq. (A3) guarantees no carries from |s j+k−1 and beyond. Therefore we can directly compute the carry from |s j+k−2 into |s j+k−1 without worrying about additional carries. This eliminates the need for an ancilla qubit |c j+k−2 , and hence the number of carry qubits needed is k − 2. To further reduce the number of ancilla qubits, we can use any sum qubit |s j = |0 during the computation. In our case, since we are processing Ω column by column, all sum qubits more significant than |s j+k−1 would be |0 . In other words, we have the last m − (j + k − 1) sum qubits usable as carry qubits in the computation described above.
As the weights are known at the time of building the circuit, the possible states that |s m can have before each addition of a state qubit |a i are also computable. Since we are adding |a i to |s m starting from the least significant bit, k equals the bit length of the maximum possible sum on |s j . . . s m after adding |a i to |s j . In other words, Therefore, the number of carry operations and additional ancilla qubits required for each addition of |a i can be determined. The term in the · in Eq. (A4) is upperbounded by u≤i, or v≤j where n max = max m j=1 i=1 nΩ i,j is the maximum number of 1's in a column of Ω. It immediately follows that the number of non-trivial carry operations (i.e. carry operations that requires M) required to add |a i to |s j . . . s m is upper-bounded by k − 2 < log 2 n max ≤ log 2 n , and the number of ancilla qubits required for the entire implementation of S is at most the upper bound for k−2, since we may use some sum qubits as carries. In other words, the number of ancilla qubits required for S grows at most logarithmically with the number of state qubits n. |a 1 + b 1 → |s 1 s 2 |a 2 + b 2 → |s 2 s 3 |a 3 + b 3 → |s 3 s 4 a 1 : • a 2 : • • a 3 : other. This reduces the CNOT gate count for QA |0 3 to 18 and the resulting circuit is reported in Fig. 19.   19. The optimized circuit for QA |0 3 used for the experiments on real quantum hardware. It requires 18 CNOT gates and 33 single qubit gates. The initial spot price is assumed to be equal to 2. The dashed boxes indicate which parts are used for A, A † , S ψ 0 , and S0. Note that due to the circuit optimization, some boxes are slightly overlapping.