Singular Fourier-Padé Series Expansion of European Option Prices (1706.06709v3)

Abstract: We apply a new numerical method, the singular Fourier-Pad\'e (SFP) method invented by Driscoll and Fornberg (2001, 2011), to price European-type options in L\'evy and affine processes. The motivation behind this application is to reduce the inefficiency of current Fourier techniques when they are used to approximate piecewise continuous (non-smooth) probability density functions. When techniques such as fast Fourier transforms and Fourier series are applied to price and hedge options with non-smooth probability density functions, they cause the Gibbs phenomenon, accordingly, the techniques converge slowly for density functions with jumps in value or derivatives. This seriously adversely affects the efficiency and accuracy of these techniques. In this paper, we derive pricing formulae and their option Greeks using the SFP method to resolve the Gibbs phenomenon and restore the global spectral convergence rate. Moreover, we show that our method requires a small number of terms to yield fast error convergence, and it is able to accurately price any European-type option deep in/out of the money and with very long/short maturities. Furthermore, we conduct an error-bound analysis of the SFP method in option pricing. This new method performs favourably in numerical experiments compared with existing techniques.

• The paper introduces the Singular Fourier-Padé method to overcome limitations of traditional Fourier techniques by addressing the Gibbs phenomenon in pricing European options.
• It details a methodology that incorporates logarithmic branch singularities into the Fourier-Padé framework, enhancing the accuracy of non-smooth PDF approximations.
• Numerical experiments across models including Black-Scholes, Variance Gamma, and Heston validate the method's fast convergence and high accuracy in option pricing.

Singular Fourier-Padé Series Expansion of European Option Prices

This paper introduces a novel numerical technique, the Singular Fourier-Padé (SFP) method, for pricing European-type options within Lévy and affine process frameworks. The core motivation lies in mitigating the inefficiencies of conventional Fourier-based methods when applied to probability density functions (PDFs) exhibiting piecewise continuity, or non-smoothness.

Addressing the Gibbs Phenomenon

Traditional Fourier techniques, including Fast Fourier Transforms (FFTs) and Fourier series, encounter challenges when dealing with non-smooth PDFs. These challenges manifest as the Gibbs phenomenon, which results in slow convergence rates and inaccuracies, particularly around discontinuities in value or derivatives. The SFP method directly addresses the Gibbs phenomenon, restoring a global spectral convergence rate.

The paper highlights two critical consequences of the Gibbs phenomenon:

1. The failure of convergence at jump discontinuities.
2. Pointwise convergence elsewhere occurs at a rate of $O(N^{-1})$, where $N$ represents the length of the Fourier partial sum.

More generally, if a function $f$ and its derivatives up to order $\nu-1$ are continuous, but $f^{(\nu)}$ is discontinuous (possessing a jump of order $\nu$), the global convergence rate is $O(N^{-\nu})$.

SFP Methodology

The SFP method builds upon the Fourier-Padé technique but incorporates logarithmic branch singularity terms. This enhancement allows for a more accurate approximation of functions with jumps and singularities. The method involves approximating a piecewise analytic function using a rational function, augmented with logarithmic terms that capture the singular behavior at jump locations.

Given a piecewise analytic function $f$ defined on the interval $[-\pi, \pi)$, with $s$ jump locations in $f$ at $t = \zeta_s\in[\pi, \pi)$, $s=1,\ldots,S,$ the complex Fourier series (CFS) representation is given by

$$f(t)=\sum_{k=-\infty}^{\infty} b_k e^{ikt},\quad\hbox{with}\quad b_k=\frac{1}{2\pi}\int_{-\pi}^{\pi} f(t)e^{-ikt} dt.$$

The transformation $z=e^{it}$, which maps the interval $[\pi,\pi)$ onto the unit circle in the complex plane, transforms the Fourier series into the following Laurent series in $z$, which can be split into

$f(z) = \sum_{k=0}^\infty{\vphantom{\sum}' b_k z^k + \sum_{k=0}^\infty{\vphantom{\sum}'b_{-k} z^{-k} = f^+(z)+f^{-}(z^{-1}),$

where the prime sums indicate that the zeroth term should be halved. The Fourier-Padé approximation of $f_{\pm}$ is comprised of polynomials

$$P_N^\pm(z)=Q_M^\pm(z)f^\pm(z)+\mathcal{O}(z^{N+M+1}), \quad z\rightarrow 0.$$

The resulting approximant is then defined as

$$\frac{P^+_N(z)}{Q^+_M(z)}+\frac{P^-_N(z^{-1})}{Q^-_M(z^{-1})}.$$

To address the approximant's limited reproduction at/around the function's jump locations, every jump in value of $f$ at $t=\zeta$ can be attributed to a logarithm of the form

$\log\left(1-\frac{z}{ e^{i\zeta}}\right).$

This logarithmic singularity in $f_\pm$ is exploited to enhance the approximation process through the condition:

$$P_N^\pm(z)+\mathcal{L}=Q_M^\pm(z)f^\pm(z)+\mathcal{O}(z^{U+1}),$$

where

$$\mathcal{L}=\sum_{s=1}^S L^{\pm}_s(z)\log\left(1-\frac{z}{ e^{i\zeta_s}}\right)$$

for some polynomials $L_s, s = 1,\ldots, S$, and $U$ is determined by S and the degrees of $P_N$, $Q_M$ and $L_s$.

Option Pricing and Greeks with SFP

The paper derives closed-form pricing formulas for European-style options and their Greeks (Delta, Gamma, Vega) using the SFP method. The method transforms the option pricing problem into a form suitable for SFP approximation. This involves representing the option price as an integral involving the payoff function and the PDF of the underlying stochastic process. The SFP method is then applied to approximate this integral, leading to efficient and accurate pricing formulas.

The SFP pricing formula of a European vanilla call option is derived as

$$V(x,K,t)=e^{-r(T-t)}\mathfrak{Re}\left[\frac{ P_N^+(z)+\sum_{s=1}^S L^+_{N_s}(z)\log\left(1-z/\varepsilon_s \right)}{ Q^+_M(z)} \right],$$

with

$$\begin{array}{rclrcl} z&=&e^{i\frac{2\pi}{d-c}(-x+\log K)}&\varepsilon_s&=&e^{i\frac{2\pi}{d-c}\zeta_s}\ P_N^+(z)&=&\sum_{n=0}^N p_{n} z^{n},&Q_M^+(z)&=&\sum_{m=0}^M q_{m} z^{m}\neq 0,\ L_{N_s}^+(z)&=&\sum_{n_s=0}^{N_s}l_{n_s} z^{n_s},& s&=&1,\ldots, S.\ \end{array}$$

Implementation and Error Analysis

The paper provides a detailed algorithm for computing the polynomial coefficients required in the SFP method. It also includes an error analysis, demonstrating that the total error in option pricing can be made arbitrarily small by choosing a suitably large integration interval. The error analysis encompasses integration truncation error, characteristic function approximation error, and SFP series truncation error.

Numerical Validation

The effectiveness of the SFP method is validated through a series of numerical experiments. These experiments cover a range of stochastic models, including the Black-Scholes-Merton (BSM), Variance Gamma (VG), CGMY, and Heston models. The SFP method consistently demonstrates fast error convergence, high accuracy, and the ability to handle options with various characteristics (deep in/out of the money, long/short maturities). The numerical results also confirm the method's global spectral convergence rate, even when dealing with non-smooth PDFs.

Conclusion

The SFP method presents a robust and efficient approach to pricing European-type options, particularly in situations where traditional Fourier methods struggle due to non-smooth PDFs. The method's ability to mitigate the Gibbs phenomenon, coupled with its fast convergence and high accuracy, makes it a valuable tool for financial modeling and risk management. Future research directions include extending the SFP method to handle more complex option types, such as spread options and options with early exercise features.