Skip to content

Monte Carlo methods: Option pricing

Overview

Many financial instruments require an estimate of the average of some function of a stochastic variable within a window of time. To compute this average, one can use Monte Carlo methods to perform many simulations of the stochastic process over the time window, evaluate the function (which can potentially depend on the path taken by the stochastic variable during the entire window), and numerically estimate the average. While the setup and details of the problems may vary from one use case to another, the underlying methods are often quite similar. As an archetypal example of this problem, we will focus on the problem of pricing derivatives, such as options, but we remark that many of these results can be carried over to other use cases, such as computing Greeks, credit valuation adjustments, value at risk, etc.

Derivatives are financial instruments that, roughly speaking, allow the parties involved to benefit when an asset (such as a stock) increases or decreases in value, but without having to already hold the asset itself. One type of derivative—called an "option"—is a contract that permits the holder to either purchase ("call option") or sell ("put option") an underlying asset at a fixed, predetermined price (the "strike price") at or prior to some predetermined time in the future ("the exercise window"). The seller of the option is obligated to either sell or buy the asset, should the holder choose to exercise the option.

How, then, should one decide on a price for the option (i.e., the amount the holder must pay for the contract, not the strike price)? The well-known Black–Scholes (or Black–Scholes–Merton) model provides one approach to pricing options, making a few assumptions about the underlying assets and the rules of the contract. More complicated options can be considered that include, for example, multiple assets in the contract (e.g., basket options), multiple possible exercise windows (e.g., Bermudan or American options), etc.

Typically, options are priced by running Monte Carlo sampling on the value of the underlying asset(s) and determining the expected profit or loss from a given position, which can be translated into a price that the purchaser must pay. Options with a larger potential downside for the seller should cost a larger amount to purchase. For more information on options and Monte Carlo methods in the context of computational finance, see [1, 2].

Actual end-to-end problem(s) solved

Suppose you want to price an option based on an underlying asset. The price of the asset is a random variable \(X\) that follows a known (or assumed) stochastic process that models the market market for the underlying asset. The option has a known payoff function \(f(X)\) (e.g., the difference between the price of the asset at each time step minus the strike price over the trajectory, or zero, whichever is larger). For options that depend on more than one underlying asset or on asset prices at multiple distinct points in time, the random variable \(X\) would represent a vector of data containing all information needed to compute the payoff. Given these inputs, the end-to-end problem is to compute the an estimate of the expected payoff \(\mathbb{E}_{X}(f(X))\) that lies within a certain error tolerance \(\epsilon\) with high probability. This quantity is then used to determine a price to charge for the option.

Using the assumed stochastic model for the price of the asset, one can develop a stochastic differential equation for the average payoff of the option. In limited cases, one can compute the average payoff analytically, as in the case of the famous Black–Scholes formula for the price of European call options, for which the 1997 Nobel Prize in economics was awarded. The Black–Scholes differential equation for the price of an asset at time \(t\) is given by

\[\begin{equation} d X_t = X_t \alpha dt + X_t \sigma d W_t, \end{equation}\]

where \(X_t\) is the price of the underlying asset at time \(t\), \(\alpha\) is a parameter known as the "drift" of the asset, \(\sigma\) is the volatility (the standard deviation of the underlying returns), and \(d W_t\) is an increment of an accompanying Brownian motion \(W_t\). Using Itô's lemma, one can derive a differential equation for the payoff function of the option at time \(t\) and, in limited cases (with several assumptions), one can solve the differential equation analytically. In practice, however, different types of contract have more complex definitions and fewer assumptions and, as a consequence, the differential equation cannot be solved analytically. Quantum approaches to numerically solving the stochastic differential equation have been proposed, including finite difference methods [3], Hamiltonian simulation [4], and quantum random walks [5], etc. For more detail on quantum approaches to solving differential equations, see the differential equations section of this document. In many real-world derivative pricing use cases, the underlying differential equation becomes intractable. Thus, the most common classical method of computing the average payoff of an option is not through solving the stochastic differential equation, but rather through Monte Carlo sampling the random process \(X\) directly. To do so, one generates a large number of price trajectories over the chosen time range, and the average payoff is computed numerically. In what follows, we will focus on quantum approaches to Monte Carlo estimation, which was pioneered in [6] and subsequently applied to several problems in finance (e.g., [7, 8, 9, 10, 11]). However, we remark that other approaches to solving this problem that do not make use of Monte Carlo methods have also been proposed (e.g., [12]), and that this is an area of active research.

If a different financial instrument is desired, such as value at risk or credit valuation adjustments, the function to be computed may be quite different, but the approach is often the same: simulate the underlying stochastic evolution several times, and estimate the desired quantity numerically.

Dominant resource cost/complexity

In [7, 8], the quantum speedup of Monte Carlo estimation from [6] is applied to solve the option pricing problem. We briefly explain the method and its dominant cost. First of all, this requires discretizing the set of values the random variable \(X\) can take, which we index by the label \(x\). Let \(N\) denote the number of values and \(n = \lceil \log_2(N) \rceil\) denote the number of qubits needed to hold the state \(\ket{x}\). The first step is to load the probabilities for the future prices of the asset into the amplitudes of a quantum state, that is, the state

\[\begin{equation} \sum_x \sqrt{p_x} \ket{x} \end{equation}\]

where \(p_x\) is the probability that \(x\) is observed in the corresponding classical Monte Carlo simulation.

Second, a subroutine is employed that computes information about the payoff function into an ancilla register using coherent arithmetic. More precisely, the angle \(\theta_x\) is computed (rounded to some finite number of bits of precision), where \(\sin(\theta_x) = \sqrt{f(x)}\). (For simplicity, here we assume \(0 \leq f(x) \leq 1\) for all \(x\), but we revisit this point later.) This yields

\[\begin{equation} \sum_x \sqrt{p_x} \ket{x}\ket{\theta_x}\,. \end{equation}\]

Third, the amplitude \(\sqrt{f(x)}\) is loaded into the amplitude of an ancilla register by applying the map \(\ket{\theta}\ket{0} \mapsto \ket{\theta}(\sin(\theta) \ket{0} + \cos(\theta)\ket{1})\). This gives

\[\begin{equation} \label{eq:option_pricing_state} \left(\sum_x \sqrt{p_xf(x)} \ket{x}\ket{\theta_x}\right)\ket{0} + \left(\sum_x \sqrt{p_x(1-f(x))} \ket{x}\ket{\theta_x}\right)\ket{1}\,. \end{equation}\]

The probability of measuring the final ancilla in \(\ket{0}\) is precisely \(\mathbb{E}_X(f(X))\). Thus, the final step is to apply quantum amplitude estimation (which requires many calls to the unitary that produces the state above) to obtain an estimate to error \(\epsilon\).

If \(0 \leq f(x) \leq 1\) does not hold, the above approach needs to be modified for example by shifting and rescaling \(f\) over a sequence of intervals of increasing length, as discussed in [6, 7]. To fit the range of \(f\) into the interval \([0,1]\), we should expect the function \(f\) will need to be scaled down by a factor on the order of the standard deviation \(\sigma = \sqrt{\mathbb{E}_X(f(X)^2)-(\mathbb{E}f(x))^2}\). Thus, to achieve error \(\epsilon\), QAE must be performed to precision \(\epsilon/\sigma\) instead of \(\epsilon\).

There are three components to the algorithm that each contribute to the resource cost:

  • Loading the distribution with amplitudes \(\sqrt{p_x}\). The gate complexity of this step is roughly the same as the time complexity of classically drawing a Monte Carlo sample, although for certain distributions it could be faster (e.g. a quadratic quantum speedup can be obtained if \(p_x\) is the stationary distribution of a Markov process [13]). Alternatively, if a functional form for \(p_x\) is known, the methods of [14] could be used to approximately prepare the state. Finally, [8] proposes using a quantum Generative Adversarial Network (qGAN), a variational ansatz, which could reduce the resources but requires a training phase.
  • Coherent arithmetic to compute the rotation angle \(\theta_x\). This depends on the complexity of the function \(f\), but can generally be accomplished in comparable gate complexity as classical arithmetic, i.e. \(\mathrm{poly}(n)\). In [15], it was shown how the payoff can instead be put directly into the amplitude, without ever computing \(\theta_x\), using quantum signal processing methods [14].
  • Quantum amplitude estimation to precision \(\epsilon/\sigma\), which requires \(\widetilde{\mathcal{O}}\left( \sigma/\epsilon \right)\) repetitions of the above two costs to achieve an \(\epsilon\)-estimate on the quantity \(\mathbb{E}_X f(X)\).

Overall, from [6, Theorem 2.5] the complexity is

\[\begin{equation} \frac{\sigma}{\epsilon}\log^{3/2}(\sigma/\epsilon)\log(\log(\sigma/\epsilon)) \cdot \mathrm{poly}(n)\,, \end{equation}\]

with the \(\mathrm{poly}(n)\) factor generally on the same order as the time required to draw and process a single classical Monte Carlo sample.

One can also consider approaches to solving this problem that rely on quantum techniques for solving differential equations, though we refer the reader to that section for more details.

Existing error corrected resource estimates

Detailed resource estimations for benchmark option-pricing problems (known as autocallable and Target Accrual Redemption Forward, or TARF) were studied in [16]. The authors studied real-world use cases and problem sizes that are relevant to current financial institutions, but on the more challenging side for classical methods. For a basket autocallable with 3 underlying assets, 5 payment days, and a knock-in put option with 20 barrier dates, the authors found that one would need about \(8000\) logical qubits, a \(T\)-depth of \(5.4\times 10^7\), and a \(T\)-count of about \(1.2\times 10^{10}\), using the most efficient methods they studied. For a TARF with 1 underlying and 26 payment dates, one needs about \(1.2 \times 10^4\) logical qubits, a \(T\)-depth of about \(8.2\times 10^7\), and a \(T\)-count of about \(9.8\times 10^9\). A follow up analysis [15] involving a quantum signal processing approach subsequently reduced these estimates to \(4.7 \times 10^3\) logical qubits, \(4.5 \times 10^7\) \(T\)-depth, and \(2.4 \times 10^9\) \(T\)-count. For comparison, classical Monte Carlo methods are roughly estimated to require 1–10 seconds and \(4 \times 10^4\) samples to achieve the same accuracy on these examples.

Similar analyses were performed in [9] for the computation of "the Greeks", which are quantities that measure the sensitivity of a derivative to various parameters. To compute the Greeks of an option, one needs to compute the derivative of the payoff function with respect to, for example, the price of the underlying. To to do this on a quantum computer, one needs to be able to estimate both the expectation of the payoff function, and have a way of computing gradients. The authors apply several quantum methods of computing gradients in order to calculate the Greeks, in addition to the quantum approaches to Monte Carlo methods used. Using a quantum gradient method to compute Greeks of an option, the authors estimate that one would need about \(1.2 \times 10^4\) logical qubits and a \(T\)-depth of around \(10^8\).

Caveats

There are many types of options and derivatives that may not be accurately captured by these simple models. Some payoff functions are path-dependent, and hence one cannot simply use the asset value at some fixed time to compute the cost, but rather the cost depends on the trajectory the random variable takes in each Monte Carlo sample.

Moreover, classical approaches to Monte Carlo sampling often allow for massive parallelization, as each simulation of the underlying asset can be done independently. By contrast, quantum algorithms for this problem require a serial approach, as the subroutines in the quantum algorithm must be run one after another without measurement and restart if the quadratic advantage is to be realized. When the slower clock speeds found in quantum devices is also taken into account, the requirements for a quantum speedup over classical methods become more stringent, as much larger problem sizes are required to achieve practical advantage. For further reading, see [17, Sec. 2.3], for example.

It is worth noting that in certain cases the number of classical samples needed to achieve error \(\epsilon\) can be reduced from the naive \(\mathcal{O}\left( \sigma^2/\epsilon^2 \right)\), cutting into the quadratic quantum speedup. In particular, quasi–Monte Carlo methods, which sample possible trajectories of the underlying assets nonrandomly can achieve a nearly quadratic speedup compared to traditional classical Monte Carlo methods, but gain an exponential dependence on the number of underlying assets ("curse of dimensionality") see [2, Chapter 5], which limits their use. The number of samples can also potentially be reduced classically via multilevel Monte Carlo methods [18]. However, when and how these methods work is delicate and must be evaluated on a case-by-case basis.

Comparable classical complexity and challenging instance sizes

Classical approaches to option pricing comprise some of the largest computational costs incurred by financial institutions. In the traditional approach to solving the option pricing problem, Monte Carlo sampling is required to simulate the evolution of the underlying asset over the time horizon of the option, and it can be slow to converge. In particular, denote the expectation value of \(f(X)\) by \(V:=\mathbb{E}_{X}(f(X))\), and the variance of \(f(X)\) by \(\sigma^2\). Classical Monte Carlo methods computes an estimate \(\hat{V}\) for \(V\) formed by averaging \(f(X)\) for \(M\) independent samples of \(X\). By Chebyshev's inequality,

\[\begin{equation} \mathrm{Pr}(|V-\hat{V}|\geq \epsilon)\leq \frac{\sigma^2}{M\epsilon^2}. \end{equation}\]

Thus, classically one needs \(M\sim\mathcal{O}(\sigma^2/\epsilon^2)\) samples to find an estimate \(\hat{V}\) within a 99% confidence interval [6].

In typical industrial scenarios, options can be priced to sufficient operational precision after roughly a few seconds of runtime, sampling as many as tens of thousands of Monte Carlo trajectories.

Alternatively, a tensor-network-based classical approach to option pricing was proposed by [19] that could lead to significant advantages over traditional classical methods in some cases.

Speedup

The classical algorithm requires \(M = \mathcal{O}\left( \sigma^2/\epsilon^2 \right)\) samples whereas the quantum algorithm requires only \(\widetilde{\mathcal{O}}\left( \sqrt{M} \right) = \widetilde{\mathcal{O}}\left( \sigma/\epsilon \right)\) samples. The gate cost of a sample is roughly the same classically and quantumly, and thus the speedup is (nearly) quadratic, inherited from the quadratic speedup of QAE.

Outlook

In [16], the authors place an upper bound on the resources required for pricing options on quantum computers, and they provide a goalpost for quantum hardware development to be able to outperform classical Monte Carlo methods. In particular, the authors estimate that a quantum device would need to be able to execute about \(10^7\) layers of \(T\)-gates per second. Moreover, the code distance for fault-tolerant implementation would need to be chosen large enough to support \(10^{10}\) total error-free logical operations. These requirements translate to a logical clock rate of about \(50\)MHz that would be needed in order to compete with current classical Monte Carlo methods. This clock speed is orders of magnitude faster than what is foreseeably possible given the current status of physical hardware and currently known methods for performing logical gates in the surface code.

While the resource requirements for pricing of derivatives are quite stringent, this is nevertheless an area of active research. For example, a new "analog" quantum representation of stochastic processes was developed in [20] that can compute \(\epsilon\)-accurate estimates of time averages (over \(T\) timesteps) of certain functions of stochastic processes in time \(O(\mathrm{polylog}(T)\epsilon^{-c})\), where \(3/2 < c <2\), an exponential speedup over classical methods in the parameter \(T\). The analog nature of their method leads to additional caveats, and finding concrete applications of this method remains an interesting open question.

Bibliography

  1. J. Hull. Options, Futures, and Other Derivatives. Pearson, 2017. ISBN 9781292212890.

  2. Paul Glasserman. Monte Carlo methods in financial engineering. Volume 53. Springer, 2004. doi:10.1007/978-0-387-21617-1.

  3. Koichi Miyamoto and Kenji Kubo. Pricing multi-asset derivatives by finite-difference method on a quantum computer. IEEE Transactions on Quantum Engineering, 3:1–25, 2021. arXiv: https://arxiv.org/abs/2109.12896. doi:10.1109/TQE.2021.3128643.

  4. Javier Gonzalez-Conde, Ángel Rodríguez-Rozas, Enrique Solano, and Mikel Sanz. Simulating option price dynamics with exponential quantum speedup. arXiv: https://arxiv.org/abs/2101.04023, 2021.

  5. Noah Linden, Ashley Montanaro, and Changpeng Shao. Quantum vs. classical algorithms for solving the heat equation. Communications in Mathematical Physics, 395(2):601–641, 2022. arXiv: https://arxiv.org/abs/2004.06516. doi:10.1007/s00220-022-04442-6.

  6. Ashley Montanaro. Quantum speedup of monte carlo methods. Proceedings of the Royal Society A, 2015. arXiv: https://arxiv.org/abs/1504.06987. doi:10.1098/rspa.2015.0301.

  7. Patrick Rebentrost, Brajesh Gupt, and Thomas R. Bromley. Quantum computational finance: monte carlo pricing of financial derivatives. Physical Review A, 98:022321, 8 2018. arXiv: https://arxiv.org/abs/1805.00109. URL: https://link.aps.org/doi/10.1103/PhysRevA.98.022321, doi:10.1103/PhysRevA.98.022321.

  8. Nikitas Stamatopoulos, Daniel J Egger, Yue Sun, Christa Zoufal, Raban Iten, Ning Shen, and Stefan Woerner. Option pricing using quantum computers. Quantum, 4:291, 2020. arXiv: https://arxiv.org/abs/1905.02666. doi:10.22331/q-2020-07-06-291.

  9. Nikitas Stamatopoulos, Guglielmo Mazzola, Stefan Woerner, and William J Zeng. Towards quantum advantage in financial market risk using quantum gradient algorithms. Quantum, 6:770, 2022. arXiv: https://arxiv.org/abs/2111.12509. doi:10.22331/q-2022-07-20-770.

  10. Jeong Yu Han and Patrick Rebentrost. Quantum advantage for multi-option portfolio pricing and valuation adjustments. arXiv: https://arxiv.org/abs/2203.04924, 2022.

  11. Stefan Woerner and Daniel J Egger. Quantum risk analysis. npj Quantum Information, 5(1):15, 2019. arXiv: https://arxiv.org/abs/1806.06893. doi:10.1038/s41534-019-0130-6.

  12. Patrick Rebentrost, Alessandro Luongo, Samuel Bosch, and Seth Lloyd. Quantum computational finance: martingale asset pricing for incomplete markets. arXiv: https://arxiv.org/abs/2209.08867, 2022.

  13. Márió Szegedy. Quantum speed-up of markov chain based algorithms. In Proceedings of the 45th IEEE Symposium on Foundations of Computer Science (FOCS), 32–41. 2004. arXiv: https://arxiv.org/abs/quant-ph/0401053. doi:10.1109/FOCS.2004.53.

  14. Sam McArdle, András Gilyén, and Mario Berta. Quantum state preparation without coherent arithmetic. arXiv: https://arxiv.org/abs/2210.14892, 2022.

  15. Nikitas Stamatopoulos and William J. Zeng. Derivative pricing using quantum signal processing. arXiv: https://arxiv.org/abs/2307.14310, 2023.

  16. Shouvanik Chakrabarti, Rajiv Krishnakumar, Guglielmo Mazzola, Nikitas Stamatopoulos, Stefan Woerner, and William J Zeng. A threshold for quantum advantage in derivative pricing. Quantum, 5:463, 2021. arXiv: https://arxiv.org/abs/2012.03819. doi:10.22331/q-2021-06-01-463.

  17. Adam Bouland, Wim van Dam, Hamed Joorati, Iordanis Kerenidis, and Anupam Prakash. Prospects and challenges of quantum finance. arXiv: https://arxiv.org/abs/2011.06492, 2020.

  18. Michael B. Giles. Multilevel monte carlo methods. Acta Numerica, 24:259–328, 2015. arXiv: https://arxiv.org/abs/1304.5472. doi:10.1017/S096249291500001X.

  19. Michael Kastoryano and Nicola Pancotti. A highly efficient tensor network algorithm for multi-asset fourier options pricing. arXiv: https://arxiv.org/abs/2203.02804, 2022.

  20. Adam Bouland, Aditi Dandapani, and Anupam Prakash. A quantum spectral method for simulating stochastic processes, with applications to monte carlo. arXiv: https://arxiv.org/abs/2303.06719, 2023.