Abstract
The Fourier-cosine expansion (COS) method is used to price European options numerically in a very efficient way. To apply the COS method, one has to specify two parameters: a truncation range for the density of the log-returns and a number of terms N to approximate the truncated density by a cosine series. How to choose the truncation range is already known. Here, we are able to find an explicit and useful bound for N as well for pricing and for the sensitivities, i.e., the Greeks Delta and Gamma, provided the density of the log-returns is smooth. We further show that the COS method has an exponential order of convergence when the density is smooth and decays exponentially. However, when the density is smooth and has heavy tails, as in the Finite Moment Log Stable model, the COS method does not have exponential order of convergence. Numerical experiments confirm the theoretical results.
Similar content being viewed by others
1 Introduction
To calibrate stock price models, it is crucial to price European options quickly because stock price models are typically calibrated to given prices of liquid call and put options by minimizing the mean-square-error between model prices and given market prices. During the optimization routine, the model prices of call and put options need to be evaluated often for different model parameters.
To compute the price of a European option, one must solve an integral involving the product of the density of the log-returns at maturity and the payoff function. However, for many financial models, the density f of the log-returns is unknown. Fortunately, the characteristic function of the log-returns is often given in closed form and can be used to obtain the density.
In their seminal paper, Fang and Oosterlee [14] proposed the COS method, which is a very efficient way to approximate the density and to compute option prices. The COS method requires two parameters: a truncation range for the density of the log-returns and a number of terms N to approximate the truncated density by a cosine series. While it is known how to choose the truncation range, see [21], the choice of N is largely based on trial and error.
The COS method has been extensively extended and applied, see [4, 15, 16, 18, 25,26,27, 34, 39, 45]. Other Fourier pricing techniques are discussed e.g., by [8, 28, 35, 36].
With respect to these papers we make the following three main contributions: we develop an explicit, useful and rigorous bound for N; we analyze the order of convergence of the COS method in detail; and we rigorously analyze how the Greeks of an option can be approximated by the COS method.
Fang and Oosterlee [14] propose to approximate the (unknown) density in three steps: (i) Truncate the density f, i.e., approximate f by a function \(f_{L}\) with finite support on some (sufficiently large) interval \([-L,L]\). (ii) Approximate \(f_{L}\) by a Fourier-cosine expansion \(\sum a_{k}e_{k}^{L}\), where \(a_{k}\) are Fourier coefficients of \(f_{L}\) and \(e_{k}^{L}\) are cosine basis functions. (iii) Approximate \(a_{k}\) by some coefficients \(c_{k}\) which can be obtained directly from the characteristic function of f. Thus, to apply the COS method, two decisions must be made: find a suitable truncation range \([-L,L]\) and identify the number N of cosine functions.
One may apply a simple triangle inequality to bound the error of the three approximations and obtain:
The first, second and third terms at the right-hand side of Inequality (1.1) correspond to approximations due to (i), (ii) and (iii), respectively.
It is well known that the series truncation error, i.e., the second term on the right-hand side of Inequality (1.1), can be bounded using integration by parts, see [7]. One contribution of this article is to use this idea in order to find an explicit and useful bound for N, provided the density of the log-returns is smooth. Our bound for N is provably large enough to ensure that the COS method converges within a predefined error tolerance. There are many financial models with smooth densities having semi-heavy tails. Examples include the Black-Scholes (BS) model, see [6], the Heston model, see [20], the Normal Inverse Gaussian (NIG) model, see [5] and the CGMY model with parameter \(Y\in (0,1)\), see [2, 3, 10, 23]. The density of the log-returns in the Variance Gamma (VG) model, see [32], is not smooth for some parameters, and our methodology cannot be applied to the VG model. We also compare the solution for finding N with another solution proposed by Aimi et al. [1].
Fang and Oosterlee [14] also analyzed the order of convergence of the COS method, focusing on the second term on the right-hand side of Inequality (1.1). They concluded that with a properly chosen truncation range, the overall error converges exponentially for smooth density functions and compares favorably to the Carr–Madan formula, see [8].
Another contribution of this article is to also consider the errors introduced by the truncation range, i.e., the errors due to (i), (ii) and (iii), and to establish upper bounds for the order of convergence of the COS method. We confirm, both theoretically and empirically, that the COS method indeed converges exponentially for smooth density functions if, in addition, the tails of the density decay at least exponentially.
However, for fat-tailed and smooth densities, such as the density of the log-returns in the Finite Moment Log Stable (FMLS) model (see [9]), the truncation error due to (i) and (iii) becomes much more relevant compared to densities with semi-heavy tails. We show theoretically that the COS method converges at least as fast as \(O(N^{-\alpha })\) for \(N\rightarrow \infty \), where \(\alpha >0\) is the Pareto tail index, e.g., for the FMLS model \(\alpha \in (1,2)\). Empirical experiments indicate that the COS method converges for such densities as fast as \(O(N^{-\alpha })\), i.e., the theoretical bound is sharp and the COS method does not converge exponentially but the order of convergence is \(\alpha \).
Greeks, also known as option sensitivities, play an important role in risk management. The Greek letters Delta or Gamma respectively represent the first and second derivatives of the price of the option with respect to the current price of the underlying asset. There are formulas in the literature on how to approximate the Delta and Gamma of the option using the COS method, see [14, 25, 38]. Another contribution of this article is to provide explicit formulas for the truncation range and the number of terms for the Greeks Delta and Gamma.
This article is structured as follows: Sect. 2 gives an overview of the technical details of the COS method. Section 3 gives explicit formulas for the truncation range and the number of terms. Section 4 analyzes the order of convergence of the COS method. In Sects. 3 and 4, we distinguish between models with semi-heavy tails and models with heavy tails. Section 5 discusses the numerical computation of the Greeks using the COS method. Section 6 contains numerical experiments that confirm the theoretical results. Section 7 concludes.
2 Overview: the COS method for option pricing
We model the stock price over time by a semimartingale \((S_{t})_{t\ge 0}\) on a filtered probability space \((\Omega ,\mathcal {F},P,(\mathcal {F}_{t})_{t\ge 0})\). The filtration \((\mathcal {F}_{t})_{t\ge 0}\) satisfies the usual conditions and \(\mathcal {F}_{0}=\{\Omega ,\emptyset \}\). We assume that there is a bank account paying continuous compound interest \(r\in \mathbb {R}\) and there is a risk-neutral measure Q. All expectations are taken under Q. All densities are risk-neutral.
There is a European option with maturity \(T>0\) and payoff \(w(S_{T})\) at T, where \(w:[0,\infty )\rightarrow \mathbb {R}\). For example, a European put option with strike \(K>0\) can be described by the payoff \(w(x)=\max (K-x,0)\), \(x\ge 0\).
In several places, we assume that the payoff function is bounded. The prices of European call options are not bounded. If we want to approximate the price of a call option with a certain error tolerance, we need only approximate the price of a put option within that error tolerance and apply the put-call parity.
Fix some \(t_{0}\in [0,T)\). The price of the European option with payoff w at time \(t_{0}\) is then given by
Since we only consider European options, we will focus on the time-0 price of the option and set \(t_{0}=0\) for the remainder of the article.
If we know the characteristic function \(\varphi _{\log (S_{T})}\) of \(\log (S_{T})\) in closed form or we are able to obtain it numerically efficiently, the COS method is able to price the European option numerically very quickly, as follows: we denote by
the centralized log-returns. The characteristic function \(\varphi \) of \(X_{T}\) is then equal to
where \(E[\log (S_{T})]=-i\varphi _{\log (S_{T})}^{\prime }(0).\) We assume that \(X_{T}\) has a density f, but the exact structure of f need not be known. Since \(E[X_T]=0\), the density of \(X_{T}\) is centered around zero and it is justified to truncate the density f on a symmetric truncation range \([-L,L]\). Define
The time-0 price of the European option with payoff w is then given by
We need some abbreviations to discuss the COS method: suppose f is \(J+1\) times continuously differentiable for \(J\ge 0\). We will approximate f by cosine functions to solve the integral at the right-hand side of Eq. (2.2) numerically. We also approximate the derivatives of f by cosine functions in order to approximate the Greeks, i.e., the sensitivities of the option, numerically.
By \(f^{(j)}\) we denote the \(j^{th}\)-derivative of f. We use the convention \(f^{(0)}\equiv f\). For \(L>0\), let \(f_{L}^{(j)}:=1_{[-L,L]}f^{(j)}\), \(j=0,\ldots ,J+1\). Suppose that \(f^{(j)}\) is integrable and vanishes at \(\pm \infty \). By integration by parts, the Fourier transform of \(f^{(j)}\) is given by
Define the basis functions
The Fourier coefficients of \(f_{L}^{(j)}\) are defined by \(a_{k}^{j}\) and approximated by \(c_{k}^{j}\), where
We also write \(a_{k}\) and \(c_{k}\) instead of \(a_{k}^{0}\) and \(c_{k}^{0}\), respectively. Intuitively, we then have
where \(\sum {}^{\prime }\) indicates that the first summand (with \(k=0\)) is weighted by one-half. A little analysis shows that
i.e., the coefficients \(c_{k}^{j}\) can be obtained explicitly if \(\varphi \) is given in closed form. Here, \(\Re (z)\) denotes the real part of a complex number z and i the imaginary unit. For \(0<M\le L\) define
To keep the notation simple, we suppress the dependence of \(a_{k}^{j}\) and \(c_{k}^{j}\) on L and the dependence of \(v_{k}\) on M. The COS method states that the time-0 price of the European option can be approximated by
The coefficients \(c_{k}\) are given in closed form when \(\varphi \) is given analytically and the coefficients \(v_{k}\) can also be computed explicitly in important cases, e.g., for plain vanilla European put or call options and digital options, see [14]. This makes the COS method numerically very efficient and robust.
In Lemma 2.1 we give the approximation in line (2.5) a precise meaning. To do so, we need a bound for the term
Integrable functions f with \(B_{f}(L)\rightarrow 0\), \(L\rightarrow \infty \) are called COS-admissible. The class of COS-admissible densities is very large; in particular, it includes bounded densities with existing first and second moments and stable densities, see [21].
Lemma 2.1
Assume \(f:\mathbb {R}\rightarrow \mathbb {R}\) is integrable and square-integrable and COS-admissible. Let \(v:\mathbb {R}\rightarrow \mathbb {R}\) be bounded, with \(|v(x)|\le K\) for all \(x\in \mathbb {R}\) and some \(K>0\). Let \(\varepsilon >0\). Let \(M>0\) so that
Define \(\xi =\sqrt{2M}K\). Let \(L\ge M\) so that
Choose N large enough so that
Then it follows that
Proof
[21, Cor. 8]. \(\square \)
Remark 2.2
Often, it is fine to choose \(M=L\), e.g., when applying the COS method to densities with semi-heavy tails. However, if the density f has heavy tails, it is usually numerically more efficient to choose L and M differently.
3 On the choice of N for smooth densities
We summarize the assumptions about the density f of the log-returns in order to find explicit expressions for M, L and N. We denote by \(C_{b}^{J+1}(\mathbb {R})\) the set of bounded functions from \(\mathbb {R}\) to \(\mathbb {R}\) which are \((J+1)\)-times, continuously differentiable with bounded derivatives. By \(\left\| .\right\| _{\infty }\) and \(\left\| .\right\| _{2}\) we denote the supremum norm and the \(L^{2}\) norm, i.e.,
Let \(\mathbb {N}_{0}=\{0\}\cup \mathbb {N}\) and \(J\in \mathbb {N}_{0}\). Let \(C_{1}>0\), \(C_{2}>0\) and \(C_{3}>0\) be suitable constants. Let \(L_{0}>0\). Assume \(f\in C_{b}^{J+1}(\mathbb {R})\). We say that f and its derivatives have semi-heavy tails if
We say that f and its derivatives have heavy tails with index \(\alpha >0\) if
We suppose that f satisfies one of the following assumptions:
Assumption A1
\(f\in C_{b}^{J+1}(\mathbb {R})\) and f and its derivatives have semi-heavy tails.
Assumption A2
\(f\in C_{b}^{J+1}(\mathbb {R})\) and f and its derivatives have heavy tails with index \(\alpha >0\).
Remark 3.1
Assume the functions \(u\mapsto |u^{j}\varphi (u)|\), \(j=0,\ldots ,J+1\), are integrable. By Fourier inversion we have that \(f(x)=\frac{1}{2\pi }\int _{\mathbb {R}}e^{iux}\varphi (u)du\). By [17, Lemma 2.8] it follows that \(f\in C_{b}^{J+1}(\mathbb {R})\).
Remark 3.2
In exceptional cases, the constants \(C_{1}\) and \(C_{2}\) are explicitly known; see Example 3.4. However, it should be pointed out that in Theorem 3.8 we obtain bounds for M, L and N for models with semi-heavy tails (e.g., BS, VG, Heston, NIG and CGMY) without knowing \(C_{1}\) or \(C_{2}\).
Remark 3.3
In Theorem 3.8 we treat models with Pareto-tails; i.e., the density for the log-returns behaves like
The right-hand side of Inequality (3.2) is obtained by differentiating the right-hand side of (3.3). We assume in Theorem 3.8 that \(C_{3}\) and \(\alpha \) are known. The exact tail-behavior of the density of the log-returns is indeed known for the stable law, in particular for the FMLS model.
Example 3.4
Let \(\sigma >0\). In the Laplace model, see [30], the centralized log-returns at maturity \(T>0\) are Laplace distributed with variance \(\sigma ^{2}T\). To ensure stock prices are finite, we need \(\sigma \sqrt{T}<\sqrt{2}\), see [19, Example 3]. It holds that
Let \(L_{0}>0\). Choose \(C_{1}=\frac{1}{\sqrt{2}\sigma \sqrt{T}}\) and \(C_{2}=\frac{\sqrt{2}}{\sigma \sqrt{T}}\). Then \(f_{\text {Lap}}^{(j)}\) satisfies Inequality (3.1).
The following lemma makes it possible to bound the series-truncation error, which depends only on the choice of N. It is known in a similar form in the literature, see e.g., Theorem 1.39 in Plonka et al. [37], Theorem 4.2 in Wright et al. [43] and Theorem 6 in Boyd [7]. It can be proven by integration by parts. It is usually stated for functions with domain \([-1,1]\) or \([0,2\pi ]\). Here, we explicitly need the dependence of the series-truncation error on the truncation range \([-L,L]\), so we give the full proof.
Lemma 3.5
Let \(J\in \mathbb {N}_{0}\). Suppose \(f\in C_{b}^{J+1}(\mathbb {R})\). It holds for \(J\ge 1\) that
and for \(J=0\) that
Proof
It holds for any \(\nu >0\) by integration by parts, see [29, Eq. (1.3)], that
For \(k\in \mathbb {N}\), we apply Equation (3.5) for \(\nu :=\frac{k\pi }{2L}\). Then it follows that
Note that \(\langle e_{k}^{L},e_{l}^{L}\rangle =L\delta _{k,l}\) for \((k,l)\ne (0,0)\) and hence
By the integral test for convergence,
which implies Eq. (3.4) for \(J\ge 1\). If \(J=0\), apply the first Equality from (3.6) and the Inequality (3.7). \(\square \)
Given \(L>0\), we need to find an upper bound for \(\left\| f^{(j)}\right\| _{\infty }\) to estimate the series truncation error by Inequality (3.4). It follows by the inverse Fourier transform and Eq. (2.3) that
Inequality (3.8) provides an explicit expression to find a bound for the term \(\Vert f^{(j)}\Vert _{\infty }\), \(j=0,\ldots ,J+1\) for several models.
Example 3.6
In the symmetric NIG model with parameters \(\alpha >0\), \(\beta =0\) and \(\delta >0\), the centralized log-returns at time \(T>0\) have density \(f_{\text {NIG}}\in C_{b}^{\infty }(\mathbb {R})\), which can be expressed in terms of the modified Bessel-function of the third kind. The characteristic function is given by
see [5] and [41, Sec. 5.3.8]. We obtain, by Inequality (3.8) and using \(\alpha ^{2}+u^{2}\ge u^{2}\), that
We need Lemma 3.7 to obtain a bound for \(B_{f}(L)\). We use the following abbreviation: the maximum of two real numbers x, y, is denoted by \(x\vee y\).
Lemma 3.7
Assume \(f\in C_{b}^{1}(\mathbb {R})\) is integrable such that \(xf^{2}(x)\rightarrow 0\), \(x\rightarrow \pm \infty \). Let \(L>0\). If
then
Proof
It holds that \(f(x)\rightarrow 0\), \(x\rightarrow \pm \infty \) and
Note that for all \(k\in \mathbb {N}\),
implying
Similarly,
Using
we arrive at
\(\square \)
In Theorem 3.8 we provide explicit formulas for M, N and L when f satisfies Assumption A1 or A2 to ensure that the COS method approximates the true price within a predefined error tolerance \(\varepsilon >0\). We also include the derivatives of f in Theorem 3.8 in order to be able to approximate the sensitivities (Greeks) of the price of the option, see Sect. 5. To approximate the time-0 price of the option by Theorem 3.8, set \(\ell =0\).
Theorem 3.8
(Find M, L and N) Let \(v:\mathbb {R}\rightarrow \mathbb {R}\) be bounded, with \(|v(x)|\le K\) for all \(x\in \mathbb {R}\) and some \(K>0\). Let \(\varepsilon >0\) be small enough. Suppose \(J\in \mathbb {N}_{0}\).
(i) Assume density f satisfies Assumption A1. For some even \(n\in \mathbb {N}\) let
where \(\mu _{n}\) is the \(n^{th}-\)moment of f, i.e., \(\mu _{n}=\frac{1}{i^{n}}\left. \frac{\partial ^{n}}{\partial u^{n}}\varphi (u)\right| _{u=0}\). Let \(\xi =\sqrt{2M}K\). If \(J\ge 1\), let \(\ell \in \{0,\ldots ,J-1\}\), \(k\in \{1,\ldots ,J-\ell \}\) and
If \(J=0\), let \(\ell =0\) and
(ii) Assume density f satisfies Assumption A2, \(J\ge 1\) and f is unimodal. Let \(M=\left( \frac{4C_{3}K}{\varepsilon \alpha }\right) ^{\frac{1}{\alpha }}\), \(\xi =\sqrt{2M}K\) and
Let \(\ell =0\) and \(k\in \{1,\ldots ,J\}\) and define N as in Equation (3.11). In both cases, i) and ii), it holds that
Proof
We start with case (i). For \(\varepsilon \) small enough, M is large enough so that Assumption A1 holds, i.e., \(L_{0}\le M=L\). It follows that
the last inequality holds true if
Further,
the last inequality holds true if
Proposition 2 in Junike and Pankrashkin [21] shows that
We then have
the last inequality holds true if
Assume \(J\ge 1\). For \(\varepsilon \) small enough, we have
because the right-hand side of Inequality (3.23) is of order \(\varepsilon ^{-\frac{1}{n}}\) while the right-hand side of Inequality (3.11) is of order \(\varepsilon ^{-(\frac{1}{n}+\frac{2}{kn}+\frac{1}{k})}\). By Inequality (3.23) it follows that
the last inequality holds true if
By Equation (3.10), M and L are of order \(\varepsilon ^{-\frac{1}{n}}\) \(.\) Hence, for \(\varepsilon \) small enough, Inequalities (3.16), (3.19), (3.22) and (3.24) are indeed satisfied because the right-hand sides of these Inequalities are of order \(\log (\varepsilon )\). By the definition of N in Inequality (3.11), we also have
By Lemma 3.5 it follows that
As \(B_{f^{(\ell )}}(L)\rightarrow 0\), \(L\rightarrow \infty \), \(f^{(\ell )}\) is COS-admissible. Inequalities (3.15), (3.18), (3.21), (3.25) and Lemma 2.1 imply Inequality (3.14). Now assume \(J=0\). By the definition of N in Inequality (3.12), Lemma 3.5 again implies Inequality (3.25). Apply Lemma 2.1 to conclude.
Next, we treat the case (ii). For \(\varepsilon \) small enough, M is large enough so that Assumption A2 holds, i.e., \(L_{0}\le M\le L\). The inequality
is satisfied by the definition of M. A little calculation shows that the definition of L implies
therefore
Since f is a unimodal density satisfying Assumption A2, f also satisfies the assumption of Lemma 3.7. To see this, note that the unimodality implies the assertion in line (3.9). Further,
Using the bound for \(B_{f}(L)\) from Lemma 3.7, we obtain
the last Inequality holds by the definition of L. For \(\varepsilon \) small enough, N is large enough and we have
Use Inequality (3.28), the definition of L and \(\frac{96}{\pi ^{2}}\le 12\), to see that
Using \(\prod _{m=1}^{j}(\alpha +m)\le (\alpha +k)^{j}\) for \(j\le k\), it follows that
the last inequality holds by Inequality (3.29). By the definition of N, we also have
By Lemma 3.5 it follows that
As \(B_{f}(L)\rightarrow 0\), \(L\rightarrow \infty \), f is COS-admissible. Apply Lemma 2.1 to conclude. \(\square \)
Remark 3.9
If v is a European put option with maturity T, K can be set to the strike of the option times \(e^{-rT}\). The error tolerance is described by \(\varepsilon \). Numerical experiments in Junike and Pankrashkin [21] suggest choosing \(n\in \{4,6,8\}\) for \(\mu _{n}\). According to Theorem 3.8, any \(k\in \{1,\ldots ,J-\ell \}\) is allowed to define N by Inequality (3.11). In the applications, one could minimize N over all admissible k. But this could be time-consuming, and in the applications, we set k to a fixed value, e.g., for the BS model, \(k=40\) is a reasonable choice, see Sect. 6.1. For other models, another choice for k might be more suitable. Bounds for \(\Vert f^{(k+1)}\Vert _{\infty }\) are explicitly known for some models, e.g., the BS, NIG and FMLS models. These bounds can also be estimated numerically, e.g., for the Heston model. Section 6 contains examples indicating that the bound for N is often sharp and very useful in applications.
Remark 3.10
If a density satisfies Assumption A1, it also satisfies Assumption A2, i.e., theoretically, case (i) in Theorem 3.8 is included in case (ii). However, in (i) we do not need to know the exact tail behavior of the density, i.e., the constants \(C_{1}\) and \(C_{2}\) from Assumption A1, in order to estimate the truncation range because we apply Markov’s inequality to find a bound for M. This approach is not applicable to densities with heavy tails because higher moments usually do not exist. In (ii), we assume the tail behavior of the density is known precisely, i.e., we have to know \(C_{3}\) and \(\alpha \) in Assumption A2 to estimate M, L or N. The constants \(C_{3}\) and \(\alpha \) are known, for example, for the FMLS model.
4 Order of convergence
Theorem 4.1 describes the order of convergence of the COS method if we allow \(N\rightarrow \infty \) and choose M and L depending on N. We describe only the asymptotic behavior of the COS method and we assume \(M=L\) in this section to keep the notation simple. We establish a bound of the order of convergence of the error of the COS method with parameters L and N, i.e.,
Let \(M=L=L(N)\). We say the error of the COS method converges with order \(\rho >0\), if there is a constant \(\kappa >0\) such that for all \(N\in \mathbb {N}\) it holds that
The error is of infinite order or converges exponentially, if Inequality (4.1) holds for any \(\rho \), see [7, Sec. 2.3]. We use big O notation: for functions \(g,h:\mathbb {N}\rightarrow \mathbb {R}\), we write \(h(N)=O(g(N))\) as \(N\rightarrow \infty \), if the absolute value of h(N) is at most a positive constant multiple of g(N) for all sufficiently large values of N.
Theorem 4.1
(Bounds for the order of convergence) Let \(v:\mathbb {R}\rightarrow \mathbb {R}\) be bounded, with \(|v(x)|\le K\) for all \(x\in \mathbb {R}\) and some \(K>0\). Assume \(J\in \mathbb {N}\).
(i) Assume density f satisfies Assumption A1. Let \(\beta \in \left( 0,\frac{J}{J+3}\right) \) and \(\gamma >0\). If \(M=L=\gamma N^{\beta }\) it holds that
(ii) Assume density f is unimodal and satisfies Assumption A2. Let \(\beta \in \left( 0,\frac{J}{J+2+2\alpha }\right) \) and \(\gamma >0\). If \(M=L=\gamma N^{\beta }\) it holds that
In both cases we have \(\text {err}(L(N),N)\rightarrow 0\), \(N\rightarrow \infty .\)
Proof
Let \(v_{L}:=1_{[-L,L]}v\). As in the proof of Corollary 8 in Junike and Pankrashkin [21], one can show that
where \(B_{f}(L)\) as in Eq. (2.6) and
We will state upper bounds for \(A_{0}\), \(A_{1}\) and \(B_{f}\) depending on the tail behaviour of f, i.e., for the different cases (i) and (ii) in Theorem 4.1. An upper bound for \(A_{2}\) can be obtained from Lemma 3.5. Note that
We now prove (i). For L large enough, we have
Further, by Inequality (3.17), it holds that
Assuming \(L\ge 1\) and applying Inequality (3.20), it follows that
Inequality (3.4) implies
where we used \(\frac{L}{N}\le \gamma \).
Let \(b_{1},\ldots ,b_{4}>0\) be suitable constants. By Inequality (4.2), it follows for N large enough that
\(\beta <\frac{J}{J+3}\) implies \(-J(1-\beta )+2\beta <0\) and the right-hand side of (4.3) converges to zero for \(N\rightarrow \infty \). We show (ii). It holds for L large enough using Inequalities (3.26) and (3.27) that
Inequality (3.4) and Assumption A2 imply for N large enough and for some suitable constants \(a_{1},\ldots ,a_{J}>0\) that
By Inequality (4.2), it holds for some suitable constants \(b_{1},\ldots ,b_{4}>0\) that
The last inequality can be seen as follows: as \(\beta \le \frac{J}{J+2+\alpha }\), it follows that
\(\square \)
Remark 4.2
The COS method converges exponentially, i.e., faster than \(O(N^{-\rho })\) as \(N\rightarrow \infty \) for any \(\rho >0\) if f is smooth and has semi-heavy tails. To see this, let \(0<\beta <1\): then,
Remark 4.3
By Theorem 4.1, the COS method converges at least like \(O(N^{-\alpha })\) as \(N\rightarrow \infty \) if f is smooth and has heavy tails with Pareto index \(\alpha >0\). Numerical experiments indicate that the COS method does not converge faster than \(O(N^{-\alpha })\) for the FMLS model, see Sect. 6.5.2, i.e., the bound in Theorem 4.1(ii) is sharp.
Remark 4.4
Theorem 4.1 cannot be applied to densities that are non-differentiable, e.g., the density of the VG model is non-differentiable for some parameters. To improve the order of convergence of the COS method if the density of the log-returns is non-differentiable, Ruijter et al. [38] apply spectral filters and consider the filter-COS method, \(\sum _{k=0}^{N}{}^{\prime }\hat{s}(\frac{k}{N})c_{k}v_{k}\), where \(\hat{s}\) is a spectral filter, i.e., a smooth function with support \([-1,1]\) and \(\hat{s}(0)=1\). For an analysis of the order of convergence and some error bounds for the filter-COS method, see [38] and references therein.
5 On the Greeks
The Greeks or sensitivities of a European option play an important role in hedging and risk management. The most important Greeks are Delta and Gamma, which are the first and second derivative of the price of a European option with respect to the current underlying price \(S_{0}>0\).
Fang and Oosterlee [14, Remark 3.2] state formulas for the approximation of Delta and Gamma by the COS method. We proof these formulas and discuss how to choose M, L and N for the Greeks.
In this section we assume \(S_{t}=S_{0}\bar{S}_{t}\) for some stochastic process \((\bar{S}_{t})_{t\ge 0}\), which does not depend on \(S_{0}\) anymore. This assumption is a very typical one, see [31, Sec. 3.1.2]. As in Sect. 2, we consider a European option with maturity \(T>0\) and payoff \(w(S_{T})\) for some \(w:[0,\infty )\rightarrow \mathbb {R}\). Let \(\eta :=E[\log (S_{T})]-\log (S_{0})\). Then \(\eta \) does not depend on \(S_{0}\). Define
The time-0 price of the European option is then given by
where, as before, f is the density of the centralized log-returns. Delta and Gamma are defined by
if the partial derivatives exist. The next lemma provides some conditions to interchange integration and differentiation in Eq. (5.2).
Lemma 5.1
Let w be bounded. Let v be defined as in Eq. (5.1). Assume \(J\in \mathbb {N}_{0}\) and density f satisfies Assumption A1. It follows that
and if \(J\ge 1\)
Proof
Let \(K>0\) such that v is bounded by K. Let \(\underline{s},\overline{s}\in \mathbb {R}\) such that \(0<\underline{s}<S_{0}<\overline{s}\). Let
Then \(x\mapsto h(x,s)\) is integrable for all \(s\in (\underline{s},\overline{s})\) and the partial derivative
exists for all \((x,s)\in \mathbb {R}\times (\underline{s},\overline{s})\). Define \(x_{0}:=L_{0}+|\log (\overline{s})|+|\log (\underline{s})|\). Then,
Let
Then \(|\frac{\partial }{\partial s}h(x,s)|\le m(x)\) for all \((x,s)\in \mathbb {R}\times (\underline{s},\overline{s})\) and \(x\mapsto m(x)\) is integrable. Interchanging differentiation and integration is allowed by the dominated convergence theorem, see e.g., [17, Lemma 2.8], and it follows that
If \(J\ge 1\), f is twice differentiable. Apply the arguments above to \(f^{(1)}\) to conclude. \(\square \)
In Theorem 5.2 we provide explicit formulas for M, N and L when f satisfies Assumption A1 to ensure that the COS method approximates the time-0 price and the Greeks Delta and Gamma within a predefined error tolerance \(\varepsilon >0\). One can use the same parameters M, N and L to obtain both the price and the Greeks. We define \(v_{k}:=\int _{-M}^{M}v(x,S_{0})e_{k}^{L}(x)dx\) as in Equation (2.4).
Theorem 5.2
(M, L and N for the time-0 price, Delta and Gamma) Let w be bounded. Let v be defined as in Equation (5.1). Let \(\varepsilon >0\) be small enough. Let \(\gamma =\min \big \{\varepsilon ,\,\varepsilon S_{0},\frac{\varepsilon S_{0}^{2}}{2}\big \}\). Suppose \(J\ge 3\). Assume density f satisfies Assumption A1. For some even \(n\in \mathbb {N}\) define
where \(\mu _{n}\) is the \(n^{th}-\)moment of f, i.e., \(\mu _{n}=\frac{1}{i^{n}}\left. \frac{\partial ^{n}}{\partial u^{n}}\varphi (u)\right| _{u=0}\). Let \(\xi =\sqrt{2M}K\), \(k\in \{1,\ldots ,J-2\}\) and
It follows that the time-0 price and the Greeks Delta and Gamma can be approximated by the COS method, i.e.,
Proof
Theorem 3.8 ensures that the time-0 price can be approximated by the COS method. By Lemma 5.1 and Theorem 3.8, it holds that
Using the triangle inequality, we can see that
\(\square \)
6 Numerical experiments
Some numerical experiments are compared with the Carr–Madan formula, see [8], which is applicable when the characteristic function of the log-returns is given in closed form and when \(E[S_{T}^{1+\gamma }]\) is finite for some \(\gamma >0\), which is the damping factor. Unless otherwise stated, we use the Carr–Madan formula with \(N=2^{17}\) terms, we set the damping factor equal to \(\gamma =0.1\) and we use a Fourier truncation range of 1200 to compute the reference prices. We implemented the Carr–Madan formula using Simpson’s rule without applying the fast Fourier transform.
All numerical experiments were performed on a modern laptop (Intel i7-10750H) using the software R and vectorized code without parallelization.
6.1 BS model
We consider the BS model with volatility \(\sigma >0\) and maturity \(T>0\). The density \(f_{\text {BS}}\) of the log-returns is normally distributed and belongs to the family of stable laws. An upper bound for \(\Vert f_{\text {BS}}^{(k+1)}\Vert _{\infty }\) can be obtained directly from Inequality (6.3) setting \(\alpha =2\) and \(c=\frac{\sigma \sqrt{T}}{\sqrt{2}}\). Let \(n\in \mathbb {N}\) be even. In the BS model, the \(n^{th}\) moment of the centralized log-returns is given by
The formula for N in Inequality (3.11) does not depend on the volatility \(\sigma \sqrt{T}\) if we set \(\ell =0\) and if we bound \(\Vert f_{\text {BS}}^{(k+1)}\Vert _{\infty }\) using Inequality (6.3). In Fig. 1, we show the dependency of N on k. At the beginning, the number of terms N decreases sharply as k increases and stabilizes approximately for \(k\ge 40\). We also see that N increases as the error tolerance \(\varepsilon \) decreases or the bound K increases.
How sharp is the bound for N in Theorems 3.8 and 5.2? We make the following experiment. Consider a put option with parameters
We set \(n=8\) and \(k=40\) to obtain \(M=L=6.94\) and \(N=218\) by Theorem 5.2. Other values for k and N are reported in Table 1. Theorem 5.2 indicates that M, L and N serves to approximate both the time-0 price and the Greeks Delta and Gamma using the COS method. This can be confirmed by an experiment: \(N_{\min }=120\) is the minimal number of terms such that the absolute differences of the approximation by the COS method and the closed form solution by the Black-Scholes formula for time-0 price, Delta and Gamma, are less than the error tolerance. N is about twice as larger as \(N_{\min }\).
How can the number of terms N be estimated if there are no closed form solutions available for the bounds of the derivatives of the density of the log-returns? We suggest solving the right-hand side of Inequality (3.8) numerically to find a bound for \(\Vert f^{(k+1)}\Vert _{\infty }\). Here, we use R’s integrate function with default values. The CPU time of the COS method and the numerical integration routine to obtain a bound for \(\Vert f^{(k+1)}\Vert _{\infty }\) are of similar magnitude; see Table 1.
The value of N does not change when using numerical integration to obtain a bound for \(\Vert f^{(k+1)}\Vert _{\infty }\) compared to the closed form solution for the bound of \(\Vert f^{(k+1)}\Vert _{\infty }\).
6.1.1 Last coefficient rule-of-thumb
By Boyd [7, Sec. 2.12], the following rule is called the Last Coefficient Rule-of-Thumb: The series truncation error approximating \(f_{L}\) by a finite cosine series is bounded by \(\sum _{k>N}|a_{k}|\). If N is large enough, the order of magnitude of the series truncation error is expected to scale like \(|a_{N}|\). [1, Sec. 3.2] propose to select the number of terms for the BS model with volatility \(\sigma >0\) and maturity T by the smallest natural number \(N_{F}\) satisfying the following inequality
where \([a,b]=[-L,L]\) is the truncation range and \(\varepsilon _{N_{F}}\) is some error tolerance. The solution to find the number of terms by (6.2) works even better than Inequality (3.11) if \(\varepsilon _{N_{F}}\) is small enough, see Table 2. However, the rule by Inequality (6.2) does not work if \(2^{-1}(b-a)\varepsilon _{N_{F}}\ge 1\) because we would then choose \(N_{F}=0\).
6.2 Heston model
There is an integral expression for the density of the log-returns in the Heston model and the density is smooth and has semi-heavy tails, see [12]. In this section, we provide numerical experiments under two different parameter sets, M1 and M2. The parameter set M1 is taken from [14] and the parameter set M2 is borrowed from [42]. We compute the time-0 prices of put options for different strikes K and different maturities T. Reference prices are obtained by the Carr–Madan formula.
Table 3 compares N obtained by Inequality (3.11) and the minimal \(N_{\min }\) such that the absolute difference of the approximation by the COS method and the reference price is less than the error tolerance \(\varepsilon \). On average, N is about six times larger than \(N_{\min }\). The CPU time using N instead of \(N_{\min }\) increases roughly by the factor three. To obtain N, we estimate \(\Vert f^{(k+1)}\Vert _{\infty }\) by solving the right-hand side of Inequality (3.8) by numeric integration using R’s function integrate with default values. The CPU times of the COS method and the numerical integration routine are of similar magnitudes.
6.3 VG model
Using the VG model as an example, this section shows that Theorem 3.8 does not help in finding the number of terms N for the COS method if the density of the log-returns is not sufficiently smooth. The density of the VG model at time \(T>0\) has semi-heavy tails, see [2, Example 7.5]. It can be expressed by means of the Whittaker function and is \((J+1)\)-times continuously differentiable if \(J+2<\frac{2T}{\nu }\), see [22].
Let \(f_{\text {VG}}\) denote the density of the log-returns in the VG model at maturity \(T>0\) with parameters \(\sigma >0\), \(\theta \in \mathbb {R}\) and \(\nu >0\). If \(T<\frac{\nu }{2}\), the density \(f_{\text {VG}}\) is unbounded and Theorem 3.8 cannot be used to find N. If \(T\in \left( \nu ,\frac{3\nu }{2}\right) \), the derivative of \(f_{\text {VG}}\) is continuous, but the second derivative of \(f_{\text {VG}}\) is unbounded, see [22, Theorem. 4.1 and Theorem. 6.1].
We can apply Theorem 3.8 with \(J=0\) if \(T\in \left( \nu ,\frac{3\nu }{2}\right) \), but the value for N is somewhat useless from a practical point of view, as the following experiment shows: Consider a European call option with the following parameters:
By Eq. (3.10), we set \(M=L=0.91\) using \(n=4\) moments. By numerically optimizing the derivative of the density \(f_{\text {VG}}\), we obtain \(\Vert f^{(1)}\Vert _{\infty }=218\), and Theorem 3.8 suggests \(N\approx 4\cdot 10^{14}\).
We calculated a reference price \(\pi _{CM}=1.809833\) using the Carr–Madan formula. Using the COS method, \(N=50\) is already sufficient to approximate the reference price within the error tolerance.
6.4 FMLS model
The stable law has been used to model stock prices since Mandelbrot [33] and Fama [13]. A representation of stable densities by special functions can be found in Zolotarev [47]. The density \(f_{\text {Stable}}\in C_{b}^{\infty }(\mathbb {R})\) of the family of stable distributions with stability parameter \(\alpha \in (0,2]\), skewness \(\beta \in [-1,1]\), scale \(c>0\) and location \(\theta \in \mathbb {R}\) has the characteristic function
see [46]. It follows by Inequality (3.8) that
where \(\Gamma \) denotes the gamma function. The density of the stable law is unimodal, see [44]. Let \(F_{\text {Stable}}\) be the cumulative distribution function for the stable density \(f_{\text {Stable}}\). By Samorodnitsky and Taqqu [40, Property 1.2.15] it holds that
where
For stable densities we therefore suggest to set \(C_{3}\) in Theorem 3.8 at least as large as \(\alpha C_{\alpha }\frac{1+|\beta |}{2}c^{\alpha }\) to obtain M and L.
The FMLS model with parameters \(\sigma >0\) and \(\alpha \in (1,2)\) describes the log-returns by a stable process with maximum negative skewness. The centralized log-returns in the FMLS model at time \(T>0\) are stably distributed with stability parameter \(\alpha \), skewness \(\beta =-1\), scale \(c=\sigma T^{\frac{1}{\alpha }}\) and location \(\theta =0\).
The density of the log-returns in the FMLS model has a heavy left tail with Pareto index \(\alpha \), i.e., the left tail decays like \(|x|^{-1-\alpha }\), \(x\rightarrow -\infty \), but the right tail decays exponentially, see [9]. This makes the FMLS very attractive from a theoretical point of view: put and call option prices and all moments of the underlying stock \(S_{T}\) exist. The expectation of the log-returns also exists, but the log-returns do not have finite variance. Fitting the FMLS model to real market data shows very satisfactory results. Carr and Wu [9] calibrated the FMLS model to real market data and estimated \(\alpha =1.5597\) and \(\sigma =0.1486\). We test the formulas for M, L and N with these values for \(\alpha \) and \(\sigma \) for a European call option with the following parameters:
The reference price is 9.7433708 by the Carr–Madan formula, and we confirm the reference price by the COS method with \(M=L=10^{5}\) and \(N=10^{7}\).
Figure 2 shows the density of the log-returns in the FMLS model and asymptotic tail behavior, i.e., the function \(x\mapsto C_{3}|x|^{-1-\alpha }\). The left tail does indeed decay like Pareto tails and the asymptotic tail behavior is very close to the density. The right tail decays faster; in fact it decays exponentially, see [9].
To apply the COS method we define by Theorem 3.8, \(M=69\), \(L=176\) and \(N=5815\) setting \(k=40\). According to Theorem 3.8, any other choice for k is possible, but another value for k does not significantly improve N.
With these parameters, the COS method prices the call option within the error tolerance in about 1.5 ms. The minimal N to stay below the error tolerance is \(N_{\min }=1200\) and the CPU time using \(N_{\min }\) is about 0.4 ms.
We also apply the Carr–Madan formula with the “default parameters”, see [8], which are also recommended by Madan and Schoutens [31, Sec. 3.1], i.e., 4096 terms, a damping factor equal to 1.5 and a Fourier truncation range of 1024. With these parameters, the Carr–Madan formula prices the option within the error tolerance in about 1.1 ms. We also implemented the Carr–Madan formula using the fast Fourier transform but we found no speed advantage, compare also with Crisóstomo [11].
The CPU time is about four times higher when using Theorem 3.8 to get N compared to the optimized value for N, which is acceptable from our point of view. The computational time of the COS method using Theorem 3.8 and that of the Carr–Madan formula with standard parameters are of similar magnitude.
6.5 Order of convergence
In this section, we empirically analyze the order of convergence of the COS method for the BS model and the FMLS model and compare the results with Theorem 4.1.
The order of convergence in (4.1) can be estimated straightforwardly by a simulation, see [24]. As
the negative slope of a straight line obtained from a log-log plot of the error \(\text {err}(L(N),N)\) against N can be used as an indicator for \(\rho \). As in Sect. 4, we assume \(M=L\) in this section.
6.5.1 BS model
In the BS model with parameter \(\sigma >0\), the density of the log-returns is arbitrarily smooth, and the tails decay even faster than exponentially. In Fig. 3 we analyze how the error of the COS method behaves in the BS model for large N and for different choices of L.
We consider a call option with parameters \(\sigma =0.2,\) \(S_{0}=K=100\), \(r=0\) and \(T=1\). We see in Fig. 3 that the COS method seems to converge exponentially at the beginning for moderate N if we choose L constant, i.e., independent of N. But for constant L, e.g., \(L=4\sigma \) or \(L=6\sigma \), the COS method does not converge for \(N\rightarrow \infty \) but the error remains constant for a certain level of N. This can be explained by Inequality (1.1): while the second term on the right-hand side of Inequality (1.1) converges to zero for \(N\rightarrow \infty \) and L fixed, the first and third terms on the right-hand side of Inequality (1.1) do not improve as N is increased but L is kept constant.
For large L, such as \(L=20\sigma \), this effect disappears somewhat because the tails of the Gaussian distribution decay so rapidly that, up to fixed-precision arithmetic,Footnote 1 the Gaussian density can be thought of as a density with finite support. Using arbitrary-precision arithmetic instead should show that even for \(L=20\sigma \), the error of the COS method does not converge to zero for \(N\rightarrow \infty \).
If we choose \(L=L(N)=\sigma \sqrt{N}\), Theorem 4.1(i) indicates that the error of the COS method converges exponentially to zero. This is confirmed empirically in Fig. 3. Other choices for L also work well, e.g., \(L=\frac{\sigma }{5}N\).
6.5.2 FMLS model
We test Theorem 4.1(ii) for the density of the log-returns in the FMLS model, which belongs to the family of stable densities and has heavy tails. For the FMLS model we use the parameters \(\alpha =1.5597\) and \(\sigma =0.1486\): these values are taken from Carr and Wu [9], who calibrated the FMLS model to real market data. We use the same reference price for the time-0 price of a European call option with maturity \(T=1\), \(S_{0}=K=100\) and \(r=0\) as in Sect. 6.4.
We compute \(L_{\text {optimal}}\) for a fixed \(N\in \mathbb {N}\) minimizing \(\text {err}(L,N)\), i.e., for each N we choose the truncation range such that the error of the COS method is minimal. In particular, for each \(N\in \{2^{m},\,m=4,5,\ldots ,20\}\) we define
where
is a sufficiently fine grid of the interval \([1,10^{6}]\).
In Fig. 4, we compute the order of convergence of the COS method for different strategies to define L. If L is constant, the COS method does not converge at all.
If we choose \(L=\frac{1}{100}N\), the empirical order of convergence is 1.57, i.e., the error behaves like \(O(N^{-1.57})\) for large N, which is very close to its theoretical bound of \(O(N^{-1.5597})\). The empirical order of convergence does not differ much if we choose \(L_{\text {optimal}}\) instead of \(L=\frac{1}{100}N\).
In particular, the numerical experiments indicate that the order of convergence is not exponential for heavy tail models even though the corresponding densities are arbitrarily smooth.
In Fig. 4 we also plot the empirical order of convergence of the Carr–Madan formula for the FMLS model with damping factor of 0.7 and a Fourier truncation range of 1200. The Figure indicates an exponential order of convergence for the Carr–Madan formula.
7 Conclusions
In this research we analyzed the COS method, which is used for efficient option pricing. The sensitivities, i.e., the Greeks, can also be efficiently approximated. The COS method requires two parameters: the truncation range \([-L,L]\) to truncate the density of the log-returns and the number of terms N to approximate the truncated density by cosine functions. We considered stock price models where the density of the log-returns is smooth and has either semi-heavy tails, i.e., the tails decay exponentially or faster, or heavy tails, i.e., Pareto tails.
In both cases, we found explicit and useful bounds for L and N and showed in numerical experiments the usefulness of these formulas in applications to obtain the time-0 price of an option and the Greeks Delta and Gamma. The densities of the log-returns are smooth for many models in finance, such as the BS, NIG, Heston and FMLS models.
If the density is not differentiable, Theorem 3.8 cannot be used to find a bound for N. If the density is only differentiable a few times, which is the case for the VG model for some parameters and short maturities, our bound for N is too large to be useful in most practical applications.
We further analyzed the order of convergence of the COS method and observed both theoretically and empirically that the models enjoy exponential convergence when the densities of the log-returns are smooth and have semi-heavy tails. However, when the density of the log-returns is smooth and has heavy tails, the error of the COS method can be bounded by \(O(N^{-\alpha })\), where \(\alpha >0\) is the Pareto tail index. This is the case, for example, for the FMLS model where \(\alpha \in (1,2)\). Numerical experiments indicate that the bound \(O(N^{-\alpha })\) is sharp.
Notes
The software package R operates with a precision of 53 bits conforming to the IEC 60559 standard, see R Core Team (2021). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/.
References
Aimi, A., Guardasoni, C., Ortiz-Gracia, L., Sanfelici, S.: Fast barrier option pricing by the COS BEM method in Heston model (with Matlab code). Comput. Methods Appl. Math. 23, 301–331 (2023)
Albin, J.M.P., Sundén, M.: On the asymptotic behaviour of Lévy processes, Part I: subexponential and exponential processes. Stoch. Process. Their Appl. 119(1), 281–304 (2009)
Asmussen, S.: On the role of skewness and kurtosis in tempered stable (CGMY) Lévy models in finance. Finance Stoch. 26(3), 383–416 (2022)
Bardgett, C., Gourier, E., Leippold, M.: Inferring volatility dynamics and risk premia from the S &P 500 and VIX markets. J. Financ. Econ. 131(3), 593–618 (2019)
Barndorff-Nielsen, O.E.: Normal inverse Gaussian distributions and stochastic volatility modelling. Scand. J. Stat. 24(1), 1–13 (1997)
Black, F., Scholes, M.: The pricing of options and corporate liabilities. J. Polit. Econ. 81(3), 637–654 (1973)
Boyd, J.P.: Chebyshev and Fourier Spectral Methods. Courier Corporation (2001)
Carr, P., Madan, D.: Option valuation using the fast Fourier transform. J. Comput. Finance 2(4), 61–73 (1999)
Carr, P., Wu, L.: The finite moment log stable process and option pricing. J. Financ. 58(2), 753–777 (2003)
Carr, P., Geman, H., Madan, D., Yor, M.: The fine structure of asset returns: an empirical investigation. J. Bus. 75(2), 305–332 (2002)
Crisóstomo, R.: Speed and biases of Fourier-based pricing choices: a numerical analysis. Int. J. Comput. Math. 95(8), 1565–1582 (2018)
Dragulescu, A.A., Yakovenko, V.M.: Probability distribution of returns in the Heston model with stochastic volatility. Quant. Finance 2(6), 443 (2002)
Fama, E.F.: The behavior of stock-market prices. J. Bus. 38(1), 34–105 (1965)
Fang, F., Oosterlee, C.W.: A novel pricing method for European options based on Fourier-cosine series expansions. SIAM J. Sci. Comput. 31(2), 826–848 (2009)
Fang, F., Oosterlee, C.W.: Pricing early-exercise and discrete barrier options by Fourier-cosine series expansions. Numer. Math. 114(1), 27 (2009)
Fang, F., Oosterlee, C.W.: A Fourier-based valuation method for Bermudan and barrier options under Heston’s model. SIAM J. Financ. Math. 2(1), 439–463 (2011)
Grubb, G.: Distributions and Operators, vol. 252. Springer, New York (2008)
Grzelak, L.A., Oosterlee, C.W.: On the Heston model with stochastic interest rates. SIAM J. Financ. Math. 2(1), 255–286 (2011)
Guillaume, F., Junike, G., Leoni, P., Schoutens, W.: Implied liquidity risk premia in option markets. Ann. Finance 15, 233–246 (2019)
Heston, S.L.: A closed-form solution for options with stochastic volatility with applications to bond and currency options. Rev. Financ. Stud. 6(2), 327–343 (1993)
Junike, G., Pankrashkin, K.: Precise option pricing by the COS method-How to choose the truncation range. Appl. Math. Comput. 421, 126935 (2022)
Küchler, U., Tappe, S.: On the shapes of bilateral Gamma densities. Stat. Probab. Lett. 78(15), 2478–2484 (2008)
Küchler, U., Tappe, S.: Tempered stable distributions and processes. Stoch. Process. Their Appl. 123(12), 4256–4293 (2013)
Leisen, D.P.J., Reimer, M.: Binomial models for option valuation-examining and improving convergence. Appl. Math. Finance 3(4), 319–346 (1996)
Leitao, Á., Oosterlee, C.W., Ortiz-Gracia, L., Bohte, S.M.: On the data-driven COS method. Appl. Math. Comput. 317, 68–84 (2018)
Liu, S., Borovykh, A., Grzelak, L.A., Oosterlee, C.W.: A neural network-based framework for financial model calibration. J. Math. Ind. 9, 1–28 (2019)
Liu, S., Oosterlee, C.W., Bohte, S.M.: Pricing options and computing implied volatilities using neural networks. Risks 7(1), 16 (2019)
Lord, R., Fang, F., Bervoets, F., Oosterlee, C.W.: A fast and accurate FFT-based method for pricing early-exercise options under Lévy processes. SIAM J. Sci. Comput. 30(4), 1678–1705 (2008)
Lyness, J.N.: Adjusted forms of the Fourier coefficient asymptotic expansion and applications in numerical quadrature. Math. Comput. 25(113), 87–104 (1971)
Madan, D.: Adapted hedging. Ann. Finance 12(3–4), 305–334 (2016)
Madan, D., Schoutens, W.: Applied Conic Finance. Cambridge University Press, Cambridge (2016)
Madan, D., Carr, P.P., Chang, E.C.: The variance gamma process and option pricing. Rev. Finance 2(1), 79–105 (1998)
Mandelbrot, B.B.: The variation of certain speculative prices. J. Bus. 38(1), 394–419 (1963)
Oosterlee, C.W., Grzelak, L.A.: Mathematical Modeling and Computation in Finance: With Exercises and Python and MATLAB Computer Codes. World Scientific (2019)
Ortiz-Gracia, L., Oosterlee, C.W.: Robust pricing of European options with wavelets and the characteristic function. SIAM J. Sci. Comput. 35(5), B1055–B1084 (2013)
Ortiz-Gracia, L., Oosterlee, C.W.: A highly efficient Shannon wavelet inverse Fourier technique for pricing European options. SIAM J. Sci. Comput. 38(1), B118–B143 (2016)
Plonka, G., Potts, D., Steidl, G., Tasche, M.: Numerical Fourier Analysis. Springer (2018)
Ruijter, M., Versteegh, M., Oosterlee, C.W.: On the application of spectral filters in a Fourier option pricing technique. J. Comput. Finance 19(1), 75–106 (2015)
Ruijter, M.J., Oosterlee, C.W.: Two-dimensional Fourier cosine series expansion method for pricing financial options. SIAM J. Sci. Comput. 34(5), B642–B671 (2012)
Samorodnitsky, G., Taqqu, M.S.: Stable Non-Gaussian Random Processes: Stochastic Models with Infinite Variance: Stochastic Modeling. Routledge (2017)
Schoutens, W.: Lévy Processes in Finance: Pricing Financial Derivatives. Wiley (2003)
Schoutens, W., Simons, E., Tistaert, J.: A perfect calibration! Now what? The best of Wilmott, pp. 281–304 (2003)
Wright, G.B., Javed, M., Montanelli, H., Trefethen, L.N.: Extension of Chebfun to periodic functions. SIAM J. Sci. Comput. 37(5), C554–C573 (2015)
Yamazato, M.: Unimodality of infinitely divisible distribution functions of class L. Ann. Probab. 6, 523–531 (1978)
Zhang, B., Oosterlee, C.W.: Efficient pricing of European-style Asian options under exponential Lévy processes based on Fourier cosine expansions. SIAM J. Financ. Math. 4(1), 399–426 (2013)
Zolotarev, V.M.: One-Dimensional Stable Distributions, vol. 65. American Mathematical Society (1986)
Zolotarev, V.M.: On representation of densities of stable laws by special functions. Theory Probab. Appl. 39(2), 354–362 (1995)
Funding
Open Access funding enabled and organized by Projekt DEAL.
Author information
Authors and Affiliations
Corresponding author
Additional information
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Junike, G. On the number of terms in the COS method for European option pricing. Numer. Math. 156, 533–564 (2024). https://doi.org/10.1007/s00211-024-01402-1
Received:
Revised:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1007/s00211-024-01402-1