1. Home
  2. Issues
  3. Volume 2, Issue 4 (2015)
  4. Option pricing in the model with stochas ...

Option pricing in the model with stochastic volatility driven by Ornstein–Uhlenbeck process. Simulation
Volume 2, Issue 4 (2015), pp. 355–369
Pub. online: 17 December 2015      Type: Research Article      Open accessOpen Access
Received
2 December 2015
Revised
10 December 2015
Accepted
10 December 2015
Published
17 December 2015

Abstract

We consider a discrete-time approximation of paths of an Ornstein–Uhlenbeck process as a mean for estimation of a price of European call option in the model of financial market with stochastic volatility. The Euler–Maruyama approximation scheme is implemented. We determine the estimates for the option price for predetermined sets of parameters. The rate of convergence of the price and an average volatility when discretization intervals tighten are determined. Discretization precision is analyzed for the case where the exact value of the price can be derived.

1 Introduction

We consider a discrete-time approximation for the price of European call option in the model of financial market with stochastic volatility driven by the Ornstein–Uhlenbeck process. An analytic expression for the price of the option is derived in [9]; however, the resulting formula is complicated and difficult to apply in most of available software. The discrete-time approximation is ready to be modeled even in the nonspecific software.
The problem of construction of discrete-time analogues for stochastic volatility models of financial markets is studied in a series of works including [5, 7, 2, 16, 1, 6, 18]. Various techniques are implemented, for example, multilevel Monte Carlo [5], conditional Monte Carlo [2, 18], exact simulation [2, 16], and Itô–Taylor approximations [7].
In most of the works, authors construct discrete-time approximations both for processes that describe the evolution of the price of asset and for processes driving the volatility of asset price. The model considered in this paper allows us to apply another approach: we only discretize the volatility process. The resulting discrete-time volatility process is then averaged in a special way and substituted into the option pricing formula. The option price is determined conditionally on the path of volatility process, and thus the conditional Monte Carlo approach is used. The rate of convergence of the option price calculated using the discrete-time volatility to the true option price for a given trajectory of volatility process is estimated.
Discretization of the model is naturally connected with the problem of discrete-time approximations to the solutions of stochastic differential equations. These matters are widely investigated and systematized in [8, 14, 17]. The simplest discrete-time approximation is the stochastic generalization of Euler approximation for deterministic differential equations proposed in [11], which is also referred to as the Euler–Maruyama scheme. Another suitable for implementation and effective method is the Milstein scheme [12]. Since the model under consideration is a diffusion with additive noise, both schemes coincide which is referred to below. It is worth noticing that Euler and Milstein schemes both belong to the class of Itô–Taylor approximations and have orders of convergence 0.5 and 1, respectively. For some diffusions, the approximation schemes can be enhanced to provide higher-order convergence, but this usually results in great increase in computation time.
Although exact simulation provides more precision compared to the Euler approximation, in this paper, we use the latter. This is motivated by the fact that the Euler approximation is cheaper in terms of computation time and by our desire to assess the rate of convergence of conditional option prices when the volatility is discretized using the Euler scheme.
This paper is structured as follows. We begin with the definition of the model under consideration and the discretization scheme used. In Section 3, the prices of the European call option are compared for discrete-time and continuous volatility processes to derive the estimate of strong convergence order. Section 4 provides numeric results of the simulation. In Section 5, we demonstrate the precision of discrete-time approximation for the case of deterministic volatility. Appendix A contains definitions and auxiliary results on discretization schemes and orders of their convergence mostly coming from [8].

2 The model and discrete approximation of volatility process

Let $\{\varOmega ,\mathcal{F},\mathbf{F}=\{{\mathcal{F}_{t}^{(B,Z)}},t\ge 0\},\mathbb{P}\}$ be a complete probability space with filtration generated by Wiener processes $\{B_{t},\hspace{0.1667em}Z_{t},\hspace{0.1667em}0\le t\le T\}$. We consider the model of the market where one risky asset is traded, its price evolves according to the geometric Brownian motion $\{S_{t},\hspace{0.2778em}0\le t\le T\},$ and its volatility is driven by a stochastic process. More precisely, the market is described by the pair of stochastic differential equations
(1)
\[dS_{t}=\mu S_{t}dt+\sigma (Y_{t})S_{t}dB_{t},\]
(2)
\[dY_{t}=-\alpha Y_{t}dt+kdZ_{t}.\]
We denote by $S_{0}=S$ and $Y_{0}=Y$ the deterministic initial values of the processes specified by Eqs. (1)–(2).
We further impose the following assumptions:
  • (C1) The Wiener processes B and Z are uncorrelated;
  • (C2) the volatility function $\sigma :\mathbb{R}\to \mathbb{R}_{+}$ is measurable, bounded away from zero by a constant c:
    \[\sigma (x)\ge c>0,\hspace{1em}x\in \mathbb{R},\]
    and satisfies the condition ${\int _{0}^{T}}{\sigma }^{2}(Y_{t})dt<\infty $ a.s.;
  • (C3) the coefficients α, μ, and k are positive.
For example, the conditions mentioned in assumption (C2) are satisfied for the measurable function $\sigma (x)$ such that $c\le {\sigma }^{2}(x)\le C$ for $0<x<T$ and some constants $0<c<C$. Moreover, given the square integrability of $\sigma (Y_{s})$, the solution of differential equation (1) is given by
(3)
\[S_{t}=S_{0}\exp \Bigg(\mu t-\frac{1}{2}{\int _{0}^{t}}{\sigma }^{2}(Y_{s})ds+{\int _{0}^{t}}\sigma (Y_{s})dB_{s}\Bigg),\]
which yields that $S_{t}$ is continuous. Hence, the product $\sigma (Y_{s})S_{t}$ is square integrable: ${\int _{0}^{T}}{\sigma }^{2}(Y_{t}){S_{t}^{2}}dt<\infty $ a.s.
The unique solution of the Langevin equation (2) $Y_{t}$ is the so called Ornstein–Uhlenbeck (OU) process. Its properties make it a suitable tool for modeling volatility in financial markets. One of the most important of the features is the mean-reversion property. The OU process is Gaussian with the following characteristics:
\[E[Y_{t}]=Y_{0}{\operatorname{e}}^{-\alpha t},\hspace{2em}\operatorname{Var}[Y_{t}]=\frac{{k}^{2}}{2\alpha }\big(1-{\operatorname{e}}^{-2\alpha t}\big).\]
Moreover, the OU process is Markov and admits the explicit representation
(4)
\[Y_{t}=Y_{0}{\operatorname{e}}^{-\alpha t}+k{\int _{0}^{t}}{\operatorname{e}}^{-\alpha (t-s)}dZ_{s}.\]
Following [9], we proceed to the risk-neutral setting characterized by the minimal martingale measure $\mathbb{Q}$. With r being the interest rate, Eqs. (1)–(2) are now in the following form (see Section 5 in [9]):
(5)
\[\begin{array}{r@{\hskip0pt}l}\displaystyle dS_{t}& \displaystyle =rS_{t}dt+\sigma (Y_{t})S_{t}d{B_{t}^{\mathbb{Q}}},\\{} \displaystyle dY_{t}& \displaystyle =-\alpha Y_{t}dt+kd{Z_{t}^{\mathbb{Q}}},\end{array}\]
where
\[\begin{array}{r@{\hskip0pt}l}\displaystyle {B_{t}^{\mathbb{Q}}}& \displaystyle =B_{t}+{\int _{0}^{t}}\frac{\mu -r}{\sigma (Y_{s})}ds\hspace{1em}\text{and}\hspace{1em}{Z_{t}^{\mathbb{Q}}}=Z_{t}\end{array}\]
are independent Wiener processes w.r.t. $\mathbb{Q}$.
This continuous-time model admits a variety of discrete-time approximations. In this paper, we apply the familiar Euler–Maruyama scheme, also referred to as the Euler scheme. The Euler–Maruyama approximation to the true solution of the Langevin equation (2) is the Markov chain ${Y}^{(m)}$ defined as follows:
  • • the partition of the interval $[0,T]$ into m equal subintervals of width $\Delta t=T/m$ is considered;
  • • the initial value of the scheme is set: ${Y_{0}^{(m)}}=Y_{0}$;
  • • ${Y_{l+1}^{(m)}}$, which we will use as a shorthand for ${Y_{(l+1)T/m}^{(m)}}$, $0\le l\le m-1$, is recursively defined by
    (6)
    \[{Y_{l+1}^{(m)}}=(1-\alpha \Delta t){Y_{l}^{(m)}}+k\Delta {Z_{l}^{\mathbb{Q}}},\]
    where $\Delta {Z_{l}^{\mathbb{Q}}}={Z_{(l+1)T/m}^{\mathbb{Q}}}-{Z_{lT/m}^{\mathbb{Q}}}$.
The continuous-time process ${Y_{t}^{(m)}}$ is a step-type process defined by
\[{Y_{t}^{(m)}}={Y_{[tm/T]T/m}^{(m)}},\hspace{1em}t\in [0,T],\]
where $[x]$ denotes an integer part of x.

3 The price of European call option

The price of European call option V in the initial time moment of in model (5) is provided by
(7)
\[V={\operatorname{e}}^{-rT}{\mathbb{E}}^{\mathbb{Q}}\big\{{\mathbb{E}}^{\mathbb{Q}}\big\{{\big({S_{T}^{\mathbb{Q}}}-K\big)}^{+}\hspace{0.1667em}\big|\hspace{0.1667em}Y_{s},0\le s\le T\big\}\big\}.\]
The inner expectation is conditional on the path of $Y_{s}$, $0\le s\le T$, and therefore, it actually is the Black–Scholes price for a model with deterministic time-dependent volatility. According to Lemma 2.1 in [13], the inner expectation in (7) has the following representation:
(8)
\[\begin{array}{r@{\hskip0pt}l}\displaystyle P& \displaystyle :={\mathbb{E}}^{\mathbb{Q}}\big\{{\big({S_{T}^{\mathbb{Q}}}-K\big)}^{+}\hspace{0.1667em}\big|\hspace{0.1667em}Y_{s},0\le s\le T\big\}={\operatorname{e}}^{\ln S+rT}\varPhi (d_{1})-K\varPhi (d_{2})\\{} & \displaystyle :={\operatorname{e}}^{\ln S+rT}\varPhi \bigg(\frac{\ln S+(r+\frac{1}{2}{\bar{\sigma }}^{2})T-\ln K}{\bar{\sigma }\sqrt{T}}\bigg)\\{} & \displaystyle \hspace{1em}-K\varPhi \bigg(\frac{\ln S+(r-\frac{1}{2}{\bar{\sigma }}^{2})T-\ln K}{\bar{\sigma }\sqrt{T}}\bigg),\end{array}\]
where $\bar{\sigma }:=\sqrt{\frac{1}{T}{\int _{0}^{T}}{\sigma }^{2}(Y_{s})ds}\ge 0$, Φ is the standard normal distribution function. The function $\bar{\sigma }$ may be viewed as the volatility averaged from the initial moment of time to maturity. The arguments of Φ in (8) are denoted as $d_{1}$ and $d_{2}$.
Our aim is to estimate the error arising as a result of approximation of the exact formula (7) by application of the Euler approximation to the process that drives volatility. Thus, we need to assess the expectation of R given by
(9)
\[R:=|P-\hat{P}_{m}|,\]
where m is the number of discretization points dividing the time interval $[0,T]$ into equal intervals, $\hat{P}_{m}$ denotes the price of the option in discrete setting calculated using the formula similar to (8):
(10)
\[\hat{P}_{m}={\operatorname{e}}^{\ln S+rT}\varPhi \big({d_{1}^{(m)}}\big)-K\varPhi \big({d_{2}^{(m)}}\big),\]
where
(11)
\[{d_{1}^{(m)}}=\frac{\ln S+(r+\frac{1}{2}{\bar{\sigma }_{m}^{2}})T-\ln K}{\bar{\sigma }_{m}\sqrt{T}},\]
(12)
\[{d_{2}^{(m)}}=\frac{\ln S+(r-\frac{1}{2}{\bar{\sigma }_{m}^{2}})T-\ln K}{\bar{\sigma }_{m}\sqrt{T}},\]
with 
(13)
\[\bar{\sigma }_{m}=\sqrt{\frac{1}{T}\sum \limits_{l=1}^{m}{\sigma }^{2}\big({Y_{l}^{(m)}}\big)\frac{T}{m}}=\sqrt{\frac{1}{m}\sum \limits_{l=1}^{m}{\sigma }^{2}\big({Y_{l}^{(m)}}\big)},\]
where ${Y_{l}^{(m)}}$ is defined in (6).
It is unlikely that we are able to find an exact or even approximate value for R. However, what really makes interest for investigation of the above bundle of models is the rate of convergence of the discrete setting to the continuous one. In order to assess the rate of convergence, the expression for an upper bound of R in terms of m needs to be derived.
Comparing (8) and (10), we can see that the approximation error arises solely due to the difference between $\bar{\sigma }$ and $\bar{\sigma }_{m}$. So, the first step would assessing the upper bound of expectation of absolute value of this difference w.r.t. m. After that, R might be expressed in terms of $R_{\sigma }:=\mathbb{E}|\bar{\sigma }-\bar{\sigma }_{m}|$.
Lemma 3.1.
Let ${\sigma }^{2}(x)$ satisfy the Hölder condition
(14)
\[\big|{\sigma }^{2}(x)-{\sigma }^{2}(y)\big|\le L|x-y{|}^{\gamma },\]
where $0<\gamma \le 1$, and L is some positive constant. Then $\mathbb{E}R_{\sigma }\le C{m}^{-0.5\gamma }$, where C is some positive constant.
Proof.
Since $\bar{\sigma }_{m}$ and $\bar{\sigma }$ are both square root functions, it is be more convenient to work with ${\bar{\sigma }_{m}^{2}}$ and ${\bar{\sigma }}^{2}$. To this end, we will use Hölder’s inequality:
\[\begin{array}{r@{\hskip0pt}l}\displaystyle \mathbb{E}|\bar{\sigma }_{m}-\bar{\sigma }|& \displaystyle =\mathbb{E}\Bigg|\sqrt{\frac{1}{T}{\int _{0}^{T}}{\sigma }^{2}(Y_{s})ds}-\sqrt{\frac{1}{m}\sum \limits_{i=1}^{m}{\sigma }^{2}\big({Y_{i}^{(m)}}\big)}\Bigg|\\{} & \displaystyle \le \mathbb{E}\Bigg|\sqrt{\Bigg|\frac{1}{T}{\int _{0}^{T}}{\sigma }^{2}(Y_{s})ds-\frac{1}{m}\sum \limits_{i=1}^{m}{\sigma }^{2}\big({Y_{i}^{(m)}}\big)\Bigg|}\Bigg|\\{} & \displaystyle \le {\Bigg(\mathbb{E}\Bigg|\frac{1}{T}{\int _{0}^{T}}{\sigma }^{2}(Y_{s})ds-\frac{1}{m}\sum \limits_{i=1}^{m}{\sigma }^{2}\big({Y_{i}^{(m)}}\big)\Bigg|\Bigg)}^{1/2}.\end{array}\]
Now we represent the integral as a sum of integrals over shorter intervals. Since the second summand does not depend on s, we may move it inside the integral sign, multiplying it by the inverse to the interval length:
\[\begin{array}{r@{\hskip0pt}l}\displaystyle \mathbb{E}|\bar{\sigma }_{m}-\bar{\sigma }|& \displaystyle \hspace{0.1667em}\le \hspace{0.1667em}{\Bigg(\mathbb{E}\Bigg|\sum \limits_{i=0}^{m-1}\Bigg(\frac{1}{T}{\int _{iT/m}^{(i+1)T/m}}{\sigma }^{2}(Y_{s})ds-\frac{1}{m}{\sigma }^{2}\big({Y_{i+1}^{(m)}}\big)\Bigg)\Bigg|\Bigg)}^{1/2}\\{} & \displaystyle \hspace{0.1667em}=\hspace{0.1667em}{\Bigg(\mathbb{E}\Bigg|\sum \limits_{i=0}^{m-1}\Bigg(\frac{1}{T}{\int _{iT/m}^{(i+1)T/m}}\hspace{-0.1667em}{\sigma }^{2}(Y_{s})ds\hspace{0.1667em}-\hspace{0.1667em}\frac{1}{m}\frac{m}{T}{\int _{iT/m}^{(i+1)T/m}}\hspace{-0.1667em}{\sigma }^{2}\big({Y_{i+1}^{(m)}}\big)ds\Bigg)\Bigg|\Bigg)}^{1/2}\\{} & \displaystyle \hspace{0.1667em}=\hspace{0.1667em}{\Bigg(\mathbb{E}\Bigg|\frac{1}{T}\sum \limits_{i=0}^{m-1}{\int _{iT/m}^{(i+1)T/m}}\big({\sigma }^{2}(Y_{s})-{\sigma }^{2}\big({Y_{i+1}^{(m)}}\big)\big)ds\Bigg|\Bigg)}^{1/2}.\end{array}\]
We apply the Hölder property of ${\sigma }^{2}(x)$:
\[\begin{array}{r@{\hskip0pt}l}& \displaystyle {\Bigg(\mathbb{E}\Bigg|\frac{1}{T}\sum \limits_{i=0}^{m-1}{\int _{iT/m}^{(i+1)T/m}}\big({\sigma }^{2}(Y_{s})-{\sigma }^{2}\big({Y_{i+1}^{(m)}}\big)\big)ds\Bigg|\Bigg)}^{1/2}\\{} & \displaystyle \hspace{1em}\le {\Bigg(\frac{L}{T}\mathbb{E}\Bigg(\sum \limits_{i=0}^{m-1}{\int _{iT/m}^{(i+1)T/m}}{\big|Y_{s}-{Y_{i+1}^{(m)}}\big|}^{\gamma }ds\Bigg)\Bigg)}^{1/2}\\{} & \displaystyle \hspace{1em}={\Bigg(\frac{L}{T}\sum \limits_{i=0}^{m-1}{\int _{iT/m}^{(i+1)T/m}}\mathbb{E}{\big|Y_{s}-{Y_{i+1}^{(m)}}\big|}^{\gamma }ds\Bigg)}^{1/2}.\end{array}\]
Recall that ${Y_{i}^{(m)}}$ is a shorthand for ${Y_{iT/m}^{(m)}}={Y_{s}^{(m)}}$, $s\in [iT/m,(i+1)T/m)$, and Proposition from the Appendix A yields that $\mathbb{E}|Y_{s}-{Y_{i+1}^{(m)}}|\le C_{1}{m}^{-1}$, where $C_{1}$ is some positive constant. We use Hölder’s inequality to derive that $\mathbb{E}|Y_{s}-{Y_{i+1}^{(m)}}{|}^{\gamma }\le {C_{1}^{\gamma }}{m}^{-\gamma }$ and arrive at
\[\mathbb{E}|\bar{\sigma }_{m}-\bar{\sigma }|\le {\bigg(\frac{L}{T}m\frac{T}{m}{C_{1}^{\gamma }}{m}^{-\gamma }\bigg)}^{1/2}=C{m}^{-\gamma /2}\]
for $C:=\sqrt{L{C_{1}^{\gamma }}}$, which proves the lemma.  □
The above lemma enables us to prove the main result of this section.
Theorem 3.1.
Let ${\sigma }^{2}(x)$ satisfy Hölder condition (14). Then $\mathbb{E}R\le D{m}^{-\gamma /2}$, where D is some positive constant.
Proof.
The function $\varPhi (x)$ has a continuous bounded derivative on $\mathbb{R}$; hence, we can use its Lipschitz property:
\[\begin{array}{r@{\hskip0pt}l}\displaystyle \mathbb{E}R& \displaystyle =\mathbb{E}|P-\hat{P}_{m}|\\{} & \displaystyle \le \mathbb{E}\big(S{\operatorname{e}}^{rT}\big|\varPhi (d_{1})-\varPhi \big({d_{1}^{(m)}}\big)\big|+K\big|\varPhi (d_{2})-\varPhi \big({d_{2}^{(m)}}\big)\big|\big)\\{} & \displaystyle \le \mathbb{E}\Big(\underset{x}{\sup }\big|f(x)\big|\big(S{\operatorname{e}}^{rT}\big|d_{1}-{d_{1}^{(m)}}\big|+K\big|d_{2}-{d_{2}^{(m)}}\big|\big)\Big),\end{array}\]
where $f(x)=\frac{1}{\sqrt{2\pi }}{\operatorname{e}}^{-\frac{{x}^{2}}{2}}$ is the density of the standard normal distribution. In the above representation,
(15)
\[\begin{array}{r@{\hskip0pt}l}\displaystyle \big|d_{1}-{d_{1}^{(m)}}\big|& \displaystyle =\bigg|\bigg(\frac{1}{\bar{\sigma }}-\frac{1}{\bar{\sigma }_{m}}\bigg)\frac{\ln S-\ln K+rT}{\sqrt{T}}+\frac{1}{2}\sqrt{T}(\bar{\sigma }-\bar{\sigma }_{m})\bigg|\end{array}\]
(16)
\[\begin{array}{r@{\hskip0pt}l}& \displaystyle \le |\bar{\sigma }-\bar{\sigma }_{m}|\bigg|\frac{1}{\bar{\sigma }\bar{\sigma }_{m}}\frac{\ln (S/K)+rT}{\sqrt{T}}+\frac{\sqrt{T}}{2}\bigg|\end{array}\]
(17)
\[\begin{array}{r@{\hskip0pt}l}& \displaystyle \le |\bar{\sigma }-\bar{\sigma }_{m}|\bigg|\frac{\ln (S/K)+rT}{{c}^{2}\sqrt{T}}+\frac{\sqrt{T}}{2}\bigg|,\end{array}\]
where c is a positive constant, and the last inequality is due to the assumption that $\sigma (x)$ is bounded away from zero for any $x\in \mathbb{R}$ (see assumption (C2)). Hence, using Lemma 3.1, we get
\[\mathbb{E}\big|d_{1}-{d_{1}^{(m)}}\big|\le C_{1}\mathbb{E}|\bar{\sigma }-\bar{\sigma }_{m}|\le D_{1}{m}^{-\gamma /2},\]
where $C_{1}:=|\frac{\ln (S/K)+rT}{{c}^{2}\sqrt{T}}+\frac{\sqrt{T}}{2}|$ and $D_{1}$ are positive constants.
Similarly, $\mathbb{E}|d_{2}-{d_{2}^{(m)}}|\le D_{2}{m}^{-\gamma /2}$, $D_{2}=const>0$, and we arrive at
(18)
\[\mathbb{E}R=\mathbb{E}|P-\hat{P}_{m}|\le \frac{1}{\sqrt{2\pi }}\big(D_{1}S{\operatorname{e}}^{rT}{m}^{-\gamma /2}+D_{2}K{m}^{-\gamma /2}\big)=D{m}^{-\gamma /2}\]
for a positive constant D.
The theorem is proved.  □

4 Numeric examples

Theorem 6.1 in [9] provides an analytic representation for the price of European call option for the stochastic volatility model under consideration. However, using it to calculate the price of an option is rather difficult and time-consuming. We further present the results of calculation of the price of European call option using simulation techniques.
The calculation process is performed in Matlab 7.9.0 and is structured as follows:
  • 1. The choice of discrete ranges of values of input parameters;
  • 2. The choice of the function $\sigma (Y_{s})$;
  • 3. For each combination of input parameters we generate 1000 trajectories of an Ornstein–Uhlenbeck process by splitting the time interval into subintervals of length $\Delta t=0.001$ and modeling values of the OU process at these points (that is, generating normally distributed variables with known mean and standard deviation using relationship (6)). For each trajectory, (10) is applied to calculate ${\bar{\sigma }_{m}^{2}}$ and the price of an option. The results for all trajectories are then averaged and discounted to provide the sample average of the price denoted by $\hat{\mathbb{E}}\hat{P}_{m}$. The average volatility over all trajectories and time interval is denoted by $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$.
Table 1.
${\sigma }^{2}(Y_{s})=a|Y_{s}|+b$
T k r K a b $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ $\hat{\mathbb{E}}\hat{P}_{m}$ $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ $\hat{\mathbb{E}}\hat{P}_{m}$
$\alpha =1$ $\alpha =100$
0.25 0.1 0 0.8 1 0 0.088 0.204 0.009 0.200
0.5 0.1 0 0.8 1 0 0.082 0.213 0.007 0.200
1 0.1 0 0.8 1 0 0.073 0.227 0.007 0.200
0.25 0.5 0 0.8 1 0 0.147 0.211 0.031 0.200
0.5 0.5 0 0.8 1 0 0.185 0.235 0.030 0.201
1 0.5 0 0.8 1 0 0.216 0.280 0.029 0.207
0.25 1 0 0.8 1 0 0.264 0.224 0.059 0.201
0.5 1 0 0.8 1 0 0.338 0.264 0.058 0.207
1 1 0 0.8 1 0 0.412 0.334 0.058 0.221
0.25 0.1 0.01 1 1 0.2 0.289 0.108 0.209 0.092
0.5 0.1 0.01 1 1 0.2 0.281 0.151 0.207 0.130
1 0.1 0.01 1 1 0.2 0.273 0.210 0.207 0.184
0.25 0.5 0.01 1 1 0.2 0.346 0.117 0.231 0.097
0.5 0.5 0.01 1 1 0.2 0.375 0.172 0.230 0.137
1 0.5 0.01 1 1 0.2 0.414 0.254 0.229 0.193
0.25 1 0.01 1 1 0.2 0.459 0.134 0.259 0.102
0.5 1 0.01 1 1 0.2 0.532 0.203 0.258 0.145
1 1 0.01 1 1 0.2 0.617 0.305 0.258 0.204
0.25 0.1 0.02 1.2 1 1 1.089 0.141 1.009 0.134
0.5 0.1 0.02 1.2 1 1 1.079 0.228 1.007 0.218
1 0.1 0.02 1.2 1 1 1.073 0.347 1.007 0.335
0.25 0.5 0.02 1.2 1 1 1.148 0.147 1.031 0.136
0.5 0.5 0.02 1.2 1 1 1.178 0.240 1.030 0.221
1 0.5 0.02 1.2 1 1 1.216 0.371 1.029 0.339
0.25 1 0.02 1.2 1 1 1.262 0.157 1.059 0.138
0.5 1 0.02 1.2 1 1 1.341 0.260 1.058 0.225
1 1 0.02 1.2 1 1 1.414 0.402 1.058 0.344
To begin with, let us recall the notation of input parameters along with ranges of values assigned to them in the process of simulation:
T – time to maturity, $T=0.25;\hspace{2.5pt}0.5;\hspace{2.5pt}1$;
k – volatility of OU process, $k=0.1;\hspace{2.5pt}0.5;\hspace{2.5pt}1$;
α – mean-reversion rate, $\alpha =1;\hspace{2.5pt}100$;
r – interest rate, $r=0;\hspace{2.5pt}0.01;\hspace{2.5pt}0.02$;
K – strike price, $K=0.8;\hspace{2.5pt}1;\hspace{2.5pt}1.2$;
$S_{0}$ – initial price of stock, $S_{0}=1$;
$Y_{0}$ – initial value of OU process, $Y_{0}=0.1$.
In order to produce numerical results, we choose the following options for the function $\sigma (Y_{s})$:
  • 1. ${\sigma }^{2}(Y_{s})=a|Y_{s}|+b$, where $a=\{0,1\}$, $b=\{0,0.2,1\}$ (Table 1);
  • 2. ${\sigma }^{2}(Y_{s})={\operatorname{e}}^{Y_{s}}+c$, $c=0.02$ (Table 2).
The results of simulations are split into groups by the mean-reversion rate α and function $\sigma (Y_{s})$. Meaningless and uninteresting results provided by some distinct combinations of inputs are ignored.
Mean-reversion of 1 corresponds to slow reverting models, and fast mean-reverting models are characterized by $\alpha =100$. Matters of speed of mean-reversion are addressed, for example, in [4].
Table 2.
${\sigma }^{2}(Y_{s})={\operatorname{e}}^{Y_{s}}+c$
T k r K $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ $\hat{\mathbb{E}}\hat{P}_{m}$ $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ $\hat{\mathbb{E}}\hat{P}_{m}$
$\alpha =1$ $\alpha =100$
0.25 0.1 0 0.8 1.113 0.303 1.024 0.297
0.5 0.1 0 0.8 1.103 0.372 1.022 0.363
1 0.1 0 0.8 1.088 0.465 1.021 0.456
0.25 0.5 0 0.8 1.135 0.305 1.025 0.297
0.5 0.5 0 0.8 1.131 0.374 1.023 0.363
1 0.5 0 0.8 1.119 0.468 1.022 0.456
0.25 1 0 0.8 1.184 0.307 1.027 0.297
0.5 1 0 0.8 1.212 0.380 1.025 0.363
1 1 0 0.8 1.238 0.478 1.024 0.456
0.25 0.1 0.01 1 1.112 0.209 1.024 0.201
0.5 0.1 0.01 1 1.103 0.291 1.022 0.281
1 0.1 0.01 1 1.086 0.401 1.021 0.390
0.25 0.5 0.01 1 1.121 0.209 1.025 0.201
0.5 0.5 0.01 1 1.129 0.294 1.023 0.281
1 0.5 0.01 1 1.128 0.405 1.022 0.390
0.25 1 0.01 1 1.178 0.213 1.026 0.201
0.5 1 0.01 1 1.206 0.299 1.025 0.281
1 1 0.01 1 1.216 0.412 1.023 0.390
0.25 0.1 0.02 1.2 1.110 0.143 1.024 0.135
0.5 0.1 0.02 1.2 1.103 0.231 1.022 0.220
1 0.1 0.02 1.2 1.087 0.349 1.021 0.338
0.25 0.5 0.02 1.2 1.133 0.145 1.025 0.135
0.5 0.5 0.02 1.2 1.128 0.233 1.023 0.220
1 0.5 0.02 1.2 1.115 0.352 1.021 0.338
0.25 1 0.02 1.2 1.162 0.147 1.027 0.135
0.5 1 0.02 1.2 1.201 0.239 1.025 0.220
1 1 0.02 1.2 1.255 0.367 1.023 0.338
We may observe that, under faster mean-reversion, the average volatility $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ and, consequently, the price of the option are lower, which is exactly what is expected from the model.
Tables 3 and 4 illustrate how the price of the option changes with the decrease of time step in discrete model.
Table 3.
${\sigma }^{2}(Y_{s})=|Y_{s}|+0.2$, $K=1$, $r=0.02$, $k=0.1$, $T=1$. Convergence
$\Delta t$ α $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ ${\bar{d}_{1}^{(m)}}$ ${\bar{d}_{2}^{(m)}}$ $\hat{\mathbb{E}}\hat{P}_{m}$
${10}^{-2}$ 1 0.272367 0.299043 −0.222056 0.213552
${10}^{-3}$ 1 0.271043 0.298506 −0.221338 0.213073
${10}^{-4}$ 1 0.272534 0.299123 −0.222179 0.213631
${10}^{-5}$ 1 0.271837 0.298822 −0.221753 0.213351
${10}^{-6}$ 1 0.271421 0.298667 −0.221560 0.213220
${10}^{-2}$ 100 0.208910 0.272291 −0.184776 0.189047
${10}^{-3}$ 100 0.206599 0.271267 −0.183264 0.188073
${10}^{-4}$ 100 0.206439 0.271196 −0.183159 0.188005
${10}^{-5}$ 100 0.206413 0.271184 −0.183142 0.187994
${10}^{-6}$ 100 0.206443 0.271198 −0.183162 0.188007
Table 4.
${\sigma }^{2}(Y_{s})={\operatorname{e}}^{Y_{s}}+0.2$, $K=1$, $r=0.02$, $k=0.1$, $T=1$. Convergence
$\Delta t$ α $\hat{\mathbb{E}}{\bar{\sigma }_{m}^{2}}$ ${\bar{d}_{1}^{(m)}}$ ${\bar{d}_{2}^{(m)}}$ $\hat{\mathbb{E}}\hat{P}_{m}$
${10}^{-2}$ 1 1.265414 0.556279 −0.519082 0.431865
${10}^{-3}$ 1 1.269504 0.579243 −0.543620 0.432472
${10}^{-4}$ 1 1.266274 0.584925 −0.549670 0.431990
${10}^{-5}$ 1 1.266030 0.566934 −0.530485 0.431948
${10}^{-6}$ 1 1.265635 0.576169 −0.540343 0.431892
${10}^{-2}$ 100 1.201083 0.566092 −0.529585 0.422128
${10}^{-3}$ 100 1.201047 0.566500 −0.530021 0.422123
${10}^{-4}$ 100 1.201026 0.566203 −0.529703 0.422120
${10}^{-5}$ 100 1.201036 0.566693 −0.530228 0.422121
${10}^{-6}$ 100 1.201023 0.566052 −0.529542 0.422119
In view of Section 3, it is also of certain interest to compare calculations obtained over one trajectory but under different discretization steps. We constructed 2000 trajectories with time-step size of ${10}^{-6}$: 1000 for the case $\alpha =1$ and 1000 for the case $\alpha =100$. These trajectories are considered to be “true” continuous-time trajectories of the Ornstein–Uhlenbeck process $Y_{t}$. The corresponding values of ${\bar{\sigma }_{m}^{2}}$ are considered to be “true” continuous-time values of ${\bar{\sigma }}^{2}$. The calculations were then performed for wider discretization intervals using the points of constructed trajectories. Thus, the samples of discretization errors for ${\bar{\sigma }_{m}^{2}}$ were derived. Probably, the estimate of ${\bar{\sigma }_{m}^{2}}$ is more valuable in such context since one would not usually calculate the price of an option over one trajectory. However, the estimate of volatility is usually derived from past data, which is in essence one distinct realization of the space of all possible scenarios.
Tables 5 and 6 provide characteristics of the samples of discretization errors. Errors are measured as a percentage of the “true” value.
Table 5.
${\sigma }^{2}(Y_{s})=|Y_{s}|+0.2$, $K=1$, $r=0.02$, $k=0.1$, $T=1$. Characteristics of sample of errors
${10}^{-2}$ ${10}^{-3}$ ${10}^{-4}$ ${10}^{-5}$
$\alpha =1$
Average 0.08710% 0.00834% 0.00081% 0.00008%
St. error 0.0000427 0.0000042 0.0000004 0
Median 0.0009575 0.0000834 0.000008 0.0000007
St. deviation 0.0013517 0.0001334 0.0000137 0.0000013
Excess −0.217306 −0.191189 −0.143295 −0.021156
Skewness 0.0492335 −0.002248 0.023124 0.0577173
Min −0.29706% −0.03669% −0.00303% −0.00036%
Max 0.52352% 0.04766% 0.00502% 0.00044%
Count 1000 1000 1000 1000
$\alpha =100$
Average 0.07790% 0.00742% 0.00083% 0.00007%
St. error 0.000043 0.0000044 0.0000004 0
Median 0.0008379 0.0000728 0.0000083 0.0000007
St. deviation 0.0013602 0.0001379 0.0000136 0.0000014
Excess −0.234452 −0.302723 −0.352995 −0.054568
Skewness −0.024765 0.0922374 0.0055451 0.0229423
Min −0.30504% −0.03231% −0.00323% −0.00037%
Max 0.46265% 0.04974% 0.00454% 0.00050%
Count 1000 1000 1000 1000
It can be seen from the tables that approximation results do not differ significantly for various time-steps. Even the widest investigated discretization interval provides acceptable precision for most applications.
Table 6.
${\sigma }^{2}(Y_{s})={\operatorname{e}}^{Y_{s}}+0.2$, $K=1$, $r=0.02$, $k=0.1$, $T=1$. Characteristics of sample of errors
${10}^{-2}$ ${10}^{-3}$ ${10}^{-4}$ ${10}^{-5}$
$\alpha =1$
Average 0.02496% 0.00268% 0.00026% 0.00002%
St. error 0.0000113 0.0000011 0.0000001 0.00000001
Median 0.0002559 0.0000266 0.0000027 0.0000002
St. deviation 0.0003584 0.0000354 0.0000035 0.0000003
Excess 0.1947561 0.1687356 −0.0576859 0.0700827
Skewness −0.1691937 −0.0097185 −0.1643507 −0.0522007
Min −0.09961% −0.00861% −0.00088% −0.00011%
Max 0.12871% 0.01464% 0.00126% 0.00013%
Count 1000 1000 1000 1000
$\alpha =100$
Average 0.02692% 0.00268% 0.00025% 0.00002%
St. error 0.0000118 0.0000012 0.0000001 0
Median 0.0002712 0.0000265 0.0000027 0.0000002
St. deviation 0.0003735 0.0000377 0.0000036 0.0000003
Excess 0.17242 0.070383 0.3383414 0.0853763
Skewness −0.0299531 −0.0195205 −0.1914745 −0.0371876
Min −0.09174% −0.01068% −0.00112% −0.00011%
Max 0.16291% 0.01411% 0.00139% 0.00014%
Count 1000 1000 1000 1000

5 Checking approximation precision in the case of deterministic volatility

In this section, we compare the option prices obtained for the Euler scheme (6) with the true prices of European call option for different sets of parameters for the case of deterministic time-dependent volatility.
The models with deterministic time-dependent volatility are the natural extension of the Black–Scholes model. The expression for the price of the option is the same as in the classical model except for the fact that, instead of constant volatility, it operates with average (or root mean square) volatility over the time interval to maturity (see, e.g., [10, 19]). Thus, the formula remains similar to (8) and (10).
It has been shown that deterministic volatility does not reflect the real-world stochastic dynamics correctly [3, 15], and such models have begun falling out of favor in the mid-1980s. The shift to stochastic volatility models was boosted by rapid development of computational tools.
Nevertheless, deterministic volatility is suitable for the purpose of our investigation since we can calculate the exact price of the option for the continuous time model.
In order to analyze the deterministic time-dependent volatility case, it looks natural to let the Brownian noise term in the definition of $Y_{t}$ vanish. Thus, we get
(19)
\[dY_{t}=-\alpha Y_{t}dt,\]
which is a familiar linear differential equation solved by
(20)
\[Y_{t}=Y_{0}{\operatorname{e}}^{-\alpha t}.\]
For the same transformation functions σ and sets of parameters as in the previous section, we calculate the prices of European call option in the continuous case using (8) and compare it with the prices of the same option calculated using (10)–(13) with
(21)
\[{Y_{l+1}^{(m)}}=(1-\alpha \Delta t){Y_{l}^{(m)}}.\]
We use the time step of 0.01 and only 10 simulations per combination of inputs. As before, all calculations are performed in Matlab 7.9.0.
Table 7 presents the results of calculations. Comparison of two approaches reveals that the Euler–Maruyama scheme provides a good approximation for the exact option price. In the case of fast mean-reversion, the results coincide when rounded to sixth digit.
Table 7.
Approximate option prices versus true option prices for deterministic volatility
T α r K a b $\hat{\mathbb{E}}\hat{P}_{m}$ $\mathbb{E}V$ $\hat{\mathbb{E}}\hat{P}_{m}$ $\mathbb{E}V$
${\sigma }^{2}(Y_{s})=a\| Y_{s}\| +b$ ${\sigma }^{2}(Y_{s})={\operatorname{e}}^{Y_{s}}+0.2$
0.25 1 0 0.8 1 0 0.203891 0.203888 0.316223 0.316220
0.5 1 0 0.8 1 0 0.211556 0.211549 0.390150 0.390147
1 1 0 0.8 1 0 0.223003 0.222994 0.490305 0.490302
0.25 1 0.01 1 1 0.2 0.107942 0.107935 0.224736 0.224733
0.5 1 0.01 1 1 0.2 0.150207 0.150199 0.312794 0.312791
1 1 0.01 1 1 0.2 0.206464 0.206457 0.429067 0.429064
0.25 1 0.02 1.2 1 1 0.141317 0.141313 0.159958 0.159954
0.5 1 0.02 1.2 1 1 0.227633 0.227629 0.253710 0.253706
1 1 0.02 1.2 1 1 0.345261 0.345257 0.379955 0.379952
0.25 100 0 0.8 1 0 0.200000 0.200000 0.309950 0.309950
0.5 100 0 0.8 1 0 0.200000 0.200000 0.382107 0.382106
1 100 0 0.8 1 0 0.200000 0.200000 0.481610 0.481610
0.25 100 0.01 1 1 0.2 0.091044 0.091044 0.217149 0.217149
0.5 100 0.01 1 1 0.2 0.128449 0.128449 0.303457 0.303457
1 100 0.01 1 1 0.2 0.181507 0.181507 0.419198 0.419198
0.25 100 0.02 1.2 1 1 0.133108 0.133108 0.152065 0.152065
0.5 100 0.02 1.2 1 1 0.217100 0.217100 0.243748 0.243748
1 100 0.02 1.2 1 1 0.333759 0.333759 0.369312 0.369312
Remark 5.1.
In this paper, we consider the price of the option at the initial time moment. However, all the above considerations are applicable for any valuation date t between the initial time moment and maturity. Some obvious changes need to be made, for example, the function $\bar{\sigma }_{t}:=\sqrt{\frac{1}{T-t}{\int _{t}^{T}}{\sigma }^{2}(Y_{s})ds}\ge 0$ needs to be introduced instead of $\bar{\sigma }$, and T needs to be substituted by $T-t$ in (8)–(13).

Appendix A. The Euler scheme: definitions and auxiliary results

The reader is advised to refer to [8], which provides in-depth study of numerical approximations of stochastic differential equations.
Consider the stochastic differential equation
(22)
\[dX_{t}=a(t,X_{t})dt+b(t,X_{t})dW_{t},\hspace{1em}t\in [t_{0},T],\]
and assume that there is a unique strong solution $X(t)$ with $X(t_{0})=X_{0}$. In order for this to be the case, certain assumptions need to be made about the functions a and b. Namely, refer to the following assumptions (assumptions (A1)–(A4) in [8], pp. 128–129):
  • A1) $a=a(t,x)$ and $b=b(t,x)$ are jointly ${L}^{2}$-measurable in $(t,x)\in [t_{0},T]\times \mathbb{R}$;
  • A2) the functions a and b satisfy the Lipschitz condition w.r.t. x, that is, there exists a constant $K>0$ such that
    \[\big|a(t,x)-a(t,y)\big|\le K|x-y|\]
    and
    \[\big|b(t,x)-b(t,y)\big|\le K|x-y|\]
    for all $t\in [t_{0},T]$ and $x,y\in \mathbb{R}$;
  • A3) there exists a constant $K>0$ such that
    \[{\big|a(t,x)\big|}^{2}\le K\big|1+|x{|}^{2}\big|\]
    and
    \[{\big|b(t,x)\big|}^{2}\le K\big|1+|x{|}^{2}\big|\]
    for all $t\in [t_{0},T]$ and $x,y\in \mathbb{R}$;
  • A4) $X_{t_{0}}$ is $\mathcal{F}_{t_{0}}$-measurable with $\mathbb{E}|X_{t_{0}}{|}^{2}<\infty $.
Let ${X_{t}^{(m)}}$ be a discretization scheme of the process $X_{t}$.
Definition.
(See [8].) We shall say that an approximating process ${X_{t}^{(m)}}$ converges in the strong sense with order $\gamma \in (0,\infty ]$ to the true process $X_{t}$ if there exists a finite constant K such that
\[\mathbb{E}\big(\big|X_{t}-{X_{t}^{(m)}}\big|\big)\le K{m}^{-\gamma }.\]
The same terminology will be applied to the functions of approximating processes.
Definition.
(See [8].) We shall say that a discrete time approximation scheme ${X_{t}^{(m)}}$ is strongly consistent if there exists a nonnegative function $c=c(m)$ with
\[\underset{m\to \infty }{\lim }c(m)=0\]
such that
\[\mathbb{E}\bigg({\bigg|\mathbb{E}\bigg(\frac{{X_{i+1}^{(m)}}-{X_{i}^{(m)}}}{T/m}\Big|\mathcal{F}_{iT/m}\bigg)-a\bigg(\frac{iT}{m},{X_{i}^{(m)}}\bigg)\bigg|}^{2}\bigg)\le c(m)\]
and
\[\mathbb{E}\bigg(\frac{m}{T}\bigg|{X_{i+1}^{(m)}}-{X_{i}^{(m)}}-\mathbb{E}\big({X_{i+1}^{(m)}}-{X_{i}^{(m)}}\big|\mathcal{F}_{iT/m}\big)-b\bigg(\frac{iT}{m},{X_{i}^{(m)}}\bigg)\Delta W_{i}{\bigg|}^{2}\bigg)\le c(m)\]
for all fixed values ${X_{i}^{(m)}}=y$ and $i=0,1,\dots ,m$.
Theorem.
(See [8], 9.6.2, p. 324.) Let assumptions (A1)–(A4) hold for (22). Then a strongly consistent equidistant-time discrete approximation ${X}^{(m)}$ of the process X on $[t_{0},T]$, with ${X_{t_{0}}^{(m)}}=X_{t_{0}}$, converges strongly to X.
Evidently, the Euler scheme ${Y}^{(m)}$ introduced to approximate Y in Section 2 satisfies all the above requirements and hence converges strongly. Moreover, it is a well-known fact that, in general, the convergence of the Euler approximation is of order 0.5. One can check these propositions using the estimates of the rate of convergence provided in [8] by the proof of Theorem 9.6.2 and Exercise 9.6.3.
However, our case is more specific since ${Y}^{(m)}$ approximates the diffusion process with additive noise, that is, $b(t,x)=k$ is constant.Hence, the following proposition holds.
Proposition.
${Y}^{(m)}$ is the Milstein scheme and thus converges strongly with order 1.
Really, the only difference in representation of ${Y}^{(m)}$ as the Milstein scheme compared to the Euler one is in the additional summand of the form
\[\frac{1}{2}b{b^{\prime }}\big({\big(\Delta {Z}^{\mathbb{Q}}\big)}^{2}-T/m\big),\]
which is identically zero for the constant function b. The Milstein scheme is known to converge with order 1 (see, e.g., [8], Theorem 10.6.3, p. 361).

References

[1] 
Andersen, L.: Simple and efficient simulation of the Heston stochastic volatility model. J. Comput. Finance 11(3), 1–42 (2008)
[2] 
Broadie, M., Kaya, O.: Exact simulation of stochastic volatility and other affine jump diffusion processes. Oper. Res. 54(2), 217–231 (2006). MR2222897. doi:10.1287/opre.1050.0247
[3] 
Brockman, P., Chowdhury, M.: Deterministic versus stochastic volatility: Implications for option pricing models. Appl. Financ. Econ. 7, 499–505 (1997)
[4] 
Fouque, J.-P., Papanicolaou, G., Sircar, K.R.: Derivatives in Financial Markets with Stochastic Volatility. Cambridge University Press, USA (2000). MR1768877
[5] 
Giles, M.B.: Multilevel Monte Carlo path simulation. Oper. Res. 56(3), 607–617 (2008). MR2436856. doi:10.1287/opre.1070.0496
[6] 
Hull, J., White, A.: The pricing of options on assets with stochastic volatilities. J. Finance 42(2), 281–300 (1987)
[7] 
Jourdain, B., Sbai, M.: High order discretization schemes for stochastic volatility models. arXiv:0908.1926 (2009)
[8] 
Kloeden, P., Platen, E.: Numerical Solution of Stochastic Differential Equations. Springer, New York (1992). MR1214374. doi:10.1007/978-3-662-12616-5
[9] 
Kuchuk-Iatsenko, S., Mishura, Y.: Pricing the European call option in the model with stochastic volatility driven by Ornstein–Uhlenbeck process. Exact formulas. Modern Stoch., Theory Appl. 2(3), 233–249 (2015)
[10] 
Lewis, A.: Option Valuation Under Stochastic Volatility: With Mathematica Code. Finance Press, Newport Beach, CA (2000). MR1742310
[11] 
Maruyama, G.: Continuous Markov processes and stochastic equations. Rend. Circ. Mat. Palermo 4, 48–90 (1955). MR0071666
[12] 
Milstein, G.: Approximate integration of stochastic differential equations. Theory Probab. Appl. 19(3), 557–562 (1975). MR0356225
[13] 
Mishura, Y., Rizhniak, G., Zubchenko, V.: European call option issued on a bond governed by a geometric or a fractional geometric Ornstein–Uhlenbeck process. Modern Stoch., Theory Appl. 1, 95–108 (2014). MR3314796. doi:10.15559/vmsta-2014.1.1.2
[14] 
Platen, E., N., B.-L.: Numerical Solution of Stochastic Differential Equations with Jumps in Finance. Springer, New York (2010). MR2723480. doi:10.1007/978-3-642-13694-8
[15] 
Sheikh, A.M.: Transaction data tests of S&P 100 call option pricing. J. Financ. Quant. Anal. 26, 459–475 (1991)
[16] 
van Haastrecht, A., Pelsser, A.: Efficient, almost exact simulation of the Heston stochastic volatility model. Int. J. Theor. Appl. Finance 31(1), 1–42 (2010). MR2646972. doi:10.1142/S0219024910005668
[17] 
Wang, H.: Monte Carlo Simulation with Applications to Finance. Chapman and Hall/CRC, Boca Raton, FL (2012). MR2962599
[18] 
Willard, G.: Calculating prices and sensitivities for path-independent derivative securities in multifactor models. J. Deriv. 5(1), 45–61 (1997)
[19] 
Wilmott, P., Howinson, S., Dewynne, I.: The Mathematics of Financial Derivatives. Cambridge University Press (1995). MR1357666. doi:10.1017/CBO9780511812545
Exit Reading PDF XML


Table of contents


Share


RSS