Stochastic Volatility

Since the appearance of the analysis by Black and Scholes, several people have tried to extend it and relax the assumptions on which the theory has been based. Merton[9] dropped the assumption of constant interest rates and showed that in this case, an option can be priced in terms of a bond price. In the same paper, Merton also showed how the Black-Scholes formula can be extended to cover the situation in which the volatility is a deterministic function of time. This is more realistic as the strong assumption of constant volatility is known not to be true[10]. Research has also been done assuming different processes for the evolution of stock prices by Merton[11], Cox and Ross[12] and Jones[13]. Cox and Ross[12] and Rubinstein[14] have solved the problem for the case when the volatility is a function of the underlying security price.

Empirical evidence investigating the distribution of stock returns has shown mixed results. Kon[15] finds that the observed distributions are consistent with stochastic volatility while Scott[16] shows that the hypothesis that stock returns are distributed independently over time can be rejected. Bodurtha and Courtadon[17] and Hull and White[18] also support the hypothesis of stochastic volatility. Considering these results, it seems reasonable to model volatility as another stochastic variable.

The Stochastic Process Followed by the Volatility

Several stochastic processes for the volatility have been considered by researchers. For example, Hull and White[19], Heston[20] and others have considered the process

$$dV = (a+bV)dt + \xi V^{1/2}Qdt$$ (66)

where $Q$ is white noise and $V = \sigma^2$. Baaquie[8], Hull and White and others have considered

$$dV = \mu Vdt + \xi VQdt$$ (67)

while Stein and Stein[21] consider

$$d\sigma = -\delta(\sigma - \theta) dt + kQdt$$ (68)

where $\delta$ and $\theta$ are constants representing the mean reversion strength and the mean value of the volatility respectively. All the processes above except for (4.3) follow the general form

$$dV = (\lambda+\mu V)dt + \xi V^\alpha Qdt$$ (69)

Hence, to keep our discussion as general as possible, we have chosen this stochastic process for the volatility. Note that this process has one more free parameter as compared to the others, namely $\alpha$. This allows for greater flexibility in the model. Most of the processes considered by researchers have been mean-reverting as there is some empirical evidence to show this. The extension of the theory to different processes for the volatility is straightforward.

The above, however, does not completely define the process followed by the volatility as there is still a possibility of a correlation between the white noise $W$ in the stock price process and $Q$, the white noise in the volatility process. Again, to keep the discussion as general as possible, we will assume that the correlation is $-1 < \rho < 1$.

The Merton-Garman Equation

The following derivation is due to Baaquie[25]

The process we are considering is
$$\begin{align} dS &= (\phi S + \sigma SW)dt \\ dV &= (\lambda + \mu V + \xi V^\alpha Q)dt \end{align}$$
where $\phi,\, \lambda,\, \mu$ and $\xi$ are constants, $V = \sigma^2$ and $W$ and $Q$ are white noise processes with correlation $-1 \le \rho \le 1$. Using Ito's lemma, we obtain the following expression for the process followed by a derivative $f_i$ dependent on the underlying security and the volatility of that security

$$\begin{split} df_i = &\left(\pdif{f_i}{t} + \phi S\pdif{f_i}... ...right) dt\\ = & \Theta_i dt + \Xi_i dX + \Psi_i dY \end{split}$$ (70)

where, in the last form, $dX = Wdt$ and $dY = Qdt$. We write it in this form to separate the stochastic and non-stochastic terms.

We now consider two different options, $f_1$ and $f_2$ on the same underlying security with strike prices and maturities given by $K_1,\, K_2,\, T_1$ and $T_2$ respectively. We form a portfolio

$$\Pi = f_1 + \Gamma_1 f_2 + \Gamma_2 S$$ (71)

so that

$$d\Pi = (\Theta_1 + \Gamma_1 \Theta_2 + \Gamma_2 \phi S)dt + ... ...a_1\Xi_2 + \Gamma_2 \sigma S)dX + (\Psi_1 + \Gamma_1 \Psi_2)dY$$ (72)

We have to get rid of the stochastic terms to ensure perfect hedging. Hence, we set
$$\begin{align} \Xi_1 + \Gamma_1 \Xi_2 + \Gamma_2\sigma S &= 0\\ \Psi_1 + \Gamma_1 \Psi_2 &= 0 \end{align}$$
to obtain
$$\begin{align} \Gamma_1 &= -\frac{\Psi_1}{\Psi_2} = -\frac{\partial f_1/\partial ... ..._1/\partial V}{\partial f_2/\partial V}\pdif{f_2}{S} - \pdif{f_1}{S} \end{align}$$

Since the portfolio is now risk-less, it must increase at the risk-free interest rate by the principle of no arbitrage. In other words, we must have

$$d\Pi = r\Pi dt$$ (73)

After expanding $\Pi$ and simplifying, we obtain

$$\begin{split} \beta(S, V, t, r) = \frac{1}{\partial f_1/\par... ...i^2 V^{2\alpha}}{2} \pdiftwo{f_2}{V} - rf_2 \right) \end{split}$$ (74)

It is important to note that $\beta$ is not a function of $K_1,\, K_2,\, T_1$ or $T_2$. This follows from the fact that the first expression is dependent only on $K_1$ and $T_1$ while the second is dependent only on $K_2$ and $T_2$. Hence, it is independent of all four variables. The term $\beta$ is referred to as the market price of volatility risk. This is because the higher the value of $\beta$, the more averse the investors are to take on the volatility risk. The reason this parameter is needed to price options with stochastic volatility and not for Black-Scholes pricing is that volatility is not traded in the market. Hence, it is not possible to perfectly hedge against the volatility even though it is possible to perfectly hedge against the underlying security price. Hence, investor risk preferences have to be considered when considering stochastic volatility or, in other words, risk-neutral valuation cannot be applied directly to the volatility since volatility is not directly traded in the market.

$\beta$ is difficult to estimate empirically and there is some evidence that it is non-zero[27]. To estimate this quantity, we consider the Cox, Ingersoll and Ross model where the consumption growth has constant correlation with the spot-asset return. This gives rise to a risk premium which is proportional to the volatility. We assume this model for simplicity as it has only the effect of redefining $\mu$ in the above equation. Henceforth, we shall assume that the market price of risk has been included in the Merton-Garman equation by redefining $\mu$. Therefore, the Merton-Garman equation for the process we are considering is

$$\pdif{f}{t} + rS\pdif{f}{S} + (\lambda + \mu V)\pdif{f}{V} +... ...}{\partial S \partial V} + \xi^2V^{2\alpha}\pdiftwo{f}{V} = rf$$ (75)

We introduce the variables $x = \ln S$ and $y = \ln V$ to simplify the calculations. In terms of these variables, the Merton-Garman equation is

$$\begin{split} &\pdif{f}{t} + \left(r-\frac{e^y}{2}\right)\pd... ... \\ &+ \xi^2 e^{2y(\alpha - 1)}\pdiftwo{f}{y} = rf \end{split}$$ (76)

For an option, we have $f(T) = \max(e^x-K, 0)$, $T$ being the maturity time. Hence, this is a final value problem.

The ``Straightforward'' Solution when $\boldsymbol{\rho = 0}$

When $\rho = 0$, the solution for any volatility process, stochastic or non-stochastic is straightforward. We make use of the theorem of Merton that the solution for a deterministic volatility process is the Black-Scholes price with the volatility variable replaced by the average volatility. We can consider the stochastic volatility case as a collection of a large number of deterministic volatility processes and the option price is then the average of the prices produced by each of the processes. In other words, if the volatility follows the generic process $V(t)$ (where $V$ may be stochastic), the option price will be given by

$$C = \int_0^\infty [SN(d_1(V)) - Ke^{-r\tau}N(d_2(V))]V_m(V)dV$$ (77)

where $V_m$ is the probability distribution function for the mean of the volatility (which is a delta function for a deterministic process) and $d_1(V)$ and $d_2(V)$ are the same variables as defined in the previous chapter.

We will give two simple examples to illustrate this. First, let us consider a deterministic process. We will choose the process

$$V = V_0e^{\mu t},\, 0 \le t \le T$$ (78)

In this case, the probability distribution function of the mean of the volatility is given by

$$V_m = \delta\left(V - V_0\frac{e^{\mu T}-1}{\mu T} \right)$$ (79)

giving us the Black-Scholes result with $\sigma$ replaced by $\sqrt{V_0\frac{e^{\mu T}-1}{\mu T}}$.

We will now consider a stochastic volatility process. We choose

$$dV = \xi Qdt,\, V(0) = V_0,\, 0 \le t \le T$$ (80)

where $Q$ represents white noise. The distribution of the mean of $V$ during the time interval $(0, T)$ is given by

$$V_m \sim N\left(V_0, \frac{\xi^2 T}{3}\right)$$ (81)

Hence, the option price is given by

$$C = \sqrt{\frac{3}{2\pi \xi^2 T}} \int_0^\infty [SN(d_1(V)) ... ...au}N(d_2(V))]\exp\left(\frac{3(V - V_0)^2}{2\xi^2 T}\right) dV$$ (82)

which can be computed numerically in a fraction of a second on a microcomputer.

However, it is usually extremely difficult to find the distribution of the mean of the volatility and hence, we usually perform a Monte-Carlo simulation of the stochastic process, finding the average volatility and the Black-Scholes price using that average volatility at each step. The prices can then be averaged over to give an estimate of the final solution. For most processes, a few tens of thousands of configurations suffice to give a solution accurate to about 4 significant figures which is more than sufficient for our purpose. Such a simulation generally takes only a couple of minutes on a Pentium 200.

We will derive Merton's theorem once we have laid the foundation of the theory for stochastic volatility. The theorem follows naturally and elegantly from the quantum mechanical formulation of the theory.

A Quantum Mechanical Formulation of the Problem

In this section, we use the quantum mechanical methods pioneered in Baaquie[8] for the volatility process

$$dV = (\lambda + \mu V + \xi V^\alpha Q)dt$$ (83)

All the steps here including the formulation of the option-pricing problem using path integration, using the Lagrangian and action to gain insight into the problem and performing the $x$-integration have been pioneered in the above paper. (In Baaquie[8], the Hamiltonian, Lagrangian and action were derived for the process $dV = (\mu V + \xi V Q)dt$. Hence, the following is an application of the method to a slightly more general stochastic process.)

The Hamiltonian for the Problem

If we define the Hamiltonian operator as

$$\begin{split} \hat{H}(x, y) =& - \left(r-\frac{e^y}{2}\right... ... &- \frac{\xi^2 e^{2y(\alpha - 1)}}{2}\pdiftwo{}{y} \end{split}$$ (84)

we obtain the Merton-Garman-Schrödinger equation

$$\pdif{f}{t} = (r + \hat{H}(x, y))f,\, f(x, y, T) = \max(e^x-K, 0)$$ (85)

which can be formally solved as

$$\begin{split} f(x, y, t) &= e^{-r\tau}\int_{-\infty}^\infty ... ...(x', T)\\ f(x, T) &= \max(e^x-K, 0),\, \tau = T-t \end{split}$$ (86)

While this looks deceptively simple, no analytic solution has been obtained for this equation. The special case $\alpha = \frac{1}{2}$ was solved using a series method by Hull and White[19] and using elementary probability techniques by Heston[20]. A solution for $\alpha = 1$ was obtained by Baaquie[8] using the path integral formalism of quantum mechanics.

The Lagrangian for the Problem

The central quantity whose knowledge is sufficient to solve the problem is the propagator

$$P(x, y, T \mid x', t) = \matel{x, y}{e^{-\hat{H}\tau}}{x'}$$ (87)

which can be conveniently handled in the Lagrangian formulation of quantum mechanics.

To determine a Lagrangian for the problem, we will have to discretise it. We discretise the time so that there are $N$ time steps.The time step is then $\epsilon = \frac{\tau}{N}$. The continuous variables $x(t)$ and $y(t)$ are then discretised to $x_i$ and $y_i$ where $0 \le i \le N$. The operator $\matel{x, y}{e^{-\hat{H}\tau}}{x', y'}$ can then be decomposed to

$$\begin{split} \matel{x, y}{e^{-\hat{H}\tau}}{x'} = &\int_{-\... ...ts \matel{x_1, y_1}{e^{-\hat{H}\epsilon}}{x_0, y_0} \end{split}$$ (88)

(where $x_N = x$, $y_N = y$ and $x_0 = x'$) in the limit $N \tendsto \infty$.

We see that if we can find $\matel{x, y}{e^{-\hat{H}\epsilon}}{x', y'}$, we can find the propagator and hence the option price. Therefore, let us look at this quantity more closely. Before we consider this quantity for the stochastic volatility case, let us consider the Black-Scholes (constant volatility) case as it is simpler and retains the essential features.

In the Black-Scholes case, we only have one variable $x$ (as $y$ is just a constant). We write

$$\matel{x}{e^{-\hat{H}_{BS}\epsilon}}{x'} = N_{BS}(\epsilon)e^{\hat{L}_{BS}\epsilon}$$ (89)

where $N(\epsilon)$ is a normalization constant. We already know the answer from (3.26). We see that
$$\begin{align} &N_{BS}(\epsilon) = \frac{1}{\sqrt{2\pi \epsilon}\sigma}\\ &\hat{... ...^2}\left(\frac{\delta x}{\epsilon} + r - \frac{\sigma^2}{2}\right)^2 \end{align}$$
where $\delta x = x-x'$.

For the stochastic volatility case, we have

$$\begin{split} N(\epsilon)e^{\hat{L}\epsilon} &= \matel{x, y}... ...t{H}\epsilon}}{p_x, p_y} \innprod{p_x, p_y}{x', y'} \end{split}$$ (90)

The Hamiltonian in the momentum basis is given by

$$\begin{split} \hat{H} =& \frac{e^y}{2}p_x^2 + \xi \rho e^{y(... ... e^{-y} - \mu -\frac{\delta y}{\epsilon}\right)ip_y \end{split}$$ (91)

Hence, we have

$$\begin{split} N(\epsilon) e^{\hat{L}\epsilon} =& \int_{-\inf... ... \left. - \xexpr ip_x - \yexpr ip_y \right) \right) \end{split}$$ (92)

(any reader who is unacquainted with multi-dimensional Gaussian integrals should take a look at appendix A) which can be done in a straightforward but tedious manner to get

$$N(\epsilon) = \frac{e^{y(1/2-\alpha)}}{2\pi\epsilon\xi \sqrt{1-\rho^2}}$$ (93)

and

$$\begin{split} \hat{L} = &-\frac{e^{2y(1-\alpha)}}{2\xi^2(1-\... ...r \yexpr \\ & - \frac{e^{-y}}{2(1-\rho^2)} \xexpr^2 \end{split}$$ (94)

which can be simplified to

$$\begin{split} \hat{L} = &-\frac{e^{-y}}{2(1-\rho^2)} \left( ... ...)^2\\ & - \frac{e^{2y(1-\alpha)}}{2\xi^2} \yexpr^2 \end{split}$$ (95)

This Lagrangian is difficult to deal with directly and hence we will consider the action to obtain an algorithm for the problem.

It should be emphasized that the above Lagrangian is only exactly correct in the limit $N \tendsto \infty$ and the complete Lagrangian may include terms of order $O(\epsilon)$ and greater apart from the above expression.

The Action for the Problem

The action, as we know, is defined as $S = \int Ldt$. Hence, the discretised version of the action is given by $S = \epsilon \sum_{i=1}^N L_i + O(\epsilon)$ where $L_i$ is the Lagrangian at time step $i$. Hence, the propagator can be written in terms of the action as
$$\begin{align} \matel{x, y}{e^{-\hat{H}\tau}}{x'} &= \int_{-\infty}^\infty dy' \m... ...-\hat{H}\tau}}{x', y'}\\ &= \lim_{N \tendsto \infty} \int DX DY e^S \end{align}$$
where we define
$$\begin{align*} DX &= \frac{e^{-y_N/2}}{\sqrt{2\pi\epsilon (1-\rho^2)}} \prod_{i=... ...ty \frac{dy_i e^{y_i(1-\alpha) }}{\sqrt{2\pi\epsilon} \xi} \right) \end{align*}$$
(again $x_0 = x',\, x_N = x,\, y_0 = y'$ and $y_N = y$). We note that the action is quadratic in $x$. This enables us to integrate over the stock price.

We also define

$$Q = \int DX e^{S_x} = \frac{e^{-y_N/2}}{\sqrt{2\pi\epsilon (... ...\frac{dx_i e^{-y_i/2}}{\sqrt{2\pi\epsilon (1-\rho^2)}} e^{S_x}$$ (96)

which is the integral of the action over the stock price.

We will now find $Q$. The $x$-dependent term in the Lagrangian is

$$L_x(i) = -\frac{e^{-y_i}}{2(1-\rho^2)} \left( \xexprsi - \frac{\rho e^{y_i(3/2-\alpha)}}{\xi}\yexpri \right)^2$$ (97)

Let

$$c_i = r - \frac{e^{y_i}}{2} - \frac{\rho e^{y_i(3/2-\alpha)}}{\xi} \yexpri$$ (98)

Hence,

$$S_x = -\frac{1}{2\epsilon(1-\rho^2)} \sum_{i=1}^N e^{-y_i} (x_i - x_{i-1} + \epsilon c_i)^2$$ (99)

We now change the variables to $z_i$ defined by

$$x_i = z_i - \epsilon \sum_{j=1}^i c_i$$ (100)

Then,

$$S_z = -\frac{1}{2\epsilon(1-\rho^2)} \sum_{i=1}^N e^{-y_i} (z_i-z_{i-1})^2$$ (101)

obtaining

$$Q = \int DX e^{S_x} = \frac{e^{-y_N/2}}{\sqrt{2\pi\epsilon (... ...\frac{dz_i e^{-y_i/2}}{\sqrt{2\pi\epsilon (1-\rho^2)}} e^{S_z}$$ (102)

All the $z_i$ integrations can be performed exactly by a process of induction. The exact procedure can be found in any textbook on path integration such as Kleinert[28] or Roepstorff[5]. (It is also treated in Baaquie[8]). We illustrate the method below.

The integration over $z_1$ is easily performed. We obtain

$$\begin{split} &\int_{-\infty}^\infty \frac{dz_1e^{-{y_1/2}}}... ...o^2)}\frac{1}{e^{y_1}+e^{y_2}} (z_2-z_0)^2 \right) \end{split}$$ (103)

The above integration can be repeatedly performed over all the variables $z_i$ to obtain

$$Q = \frac{e^{S_1}}{\sqrt{2\pi\epsilon (1-\rho^2) \sum_{i=1}^N e^{y_i}}}$$ (104)

where

$$\begin{split} S_1 &= -\frac{1}{2\epsilon (1-\rho^2) \sum_{i=... ...^{y_i(3/2-\alpha)}}{\xi} \yexpri\right)\right)^2\\ \end{split}$$ (105)

which, on taking the limit, $N \tendsto \infty$, becomes

$$S_1 = -\frac{1}{2(1-\rho^2) \omega} \left( x-x' + r\tau - \f... ... \frac{\rho \mu}{\xi}\eta + \frac{\rho \xi}{2} \zeta \right)^2$$ (106)

(for $\alpha \ne 3/2$ the term $e^{y(0)(3/2-\alpha)}$ arises from the fact that $\int_0^\tau dt e^{y(3/2-\alpha)} \frac{dy}{dt} = \int_{y(0)}^{y(\tau)} dy e^{y(3/2-\alpha)}$ and it will be easily seen that when $\alpha = 3/2$, that term is replaced by $\frac{\rho}{\xi} \left(y(\tau) - y(0)\right)$) where
$$\begin{align} \omega = &\int_0^\tau e^y dt = \int_0^\tau Vdt\\ \theta = &\int_0... ... = &\int_0^\tau e^{y(\alpha-1/2)} dt = \int_0^\tau V^{\alpha-1/2} dt \end{align}$$
when $V$ follows the random process (4.4). Hence, if we can find the joint probability density functions for $\omega,\, \theta,\, \eta,\, \zeta$ and $\nu = V^{3/2-\alpha}(0)$ with $V$ following the given stochastic process, we can get an analytic solution for the problem which will be given by

$$\matel{x, y}{e^{-\hat{H}\tau}}{x'} = \int_0^\infty d\omega \... ...1-\rho^2) \omega}} f(\omega,\, \theta,\, \eta,\, \zeta,\, \nu)$$ (107)

where $f$ is the joint probability distribution function. However, we were unable find the joint probability distribution function for the above quantities.

Hence, we retain the discrete solution which finally gives us

$$\matel{x, y}{e^{-\hat{H}\tau}}{x'} = \int DY \frac{e^{S_0+S_1}}{ \sqrt{2\pi\epsilon (1-\rho^2) \sum_{i=1}^N e^{y_i}}}$$ (108)

where $S_1$ is given above and
$$\begin{align} S_0 &= -\frac{\epsilon}{2\xi^2} \sum_{i=1}^N e^{2y_i(1-\alpha)}\ye... ...^{N-1} \frac{dy_i e^{y_i(1-\alpha)}}{\sqrt{2\pi \epsilon}\xi}\right) \end{align}$$

We are now in a position to derive a Monte-Carlo algorithm to calculate option prices with the volatility performing the stochastic process (4.4) with any correlation $-1 \le \rho \le 1$ with almost the same efficiency as the straightforward solution when $\rho = 0$. However, the method has a disadvantage in that it cannot handle lumpy dividends (for a continuous dividend yield $q$, we can just replace $r$ by $r-q$).

Merton's Theorem

We now consider the case $\rho = 0$ and derive the result we stated earlier in this chapter. In this case, the Hamiltonian, Lagrangian and action are given by
$$\begin{align} \hat{H} = &-\left(r - \frac{e^y}{2} \right) \pdif{}{x} - \left( \l... ...+ \epsilon \sum_{i=1}^N \left(r - \frac{e^{y_i}}{2}\right) \right)^2 \end{align}$$
The last expression is particularly interesting as $S_1$ which finally determines the option price is the same as that for the Black-Scholes case with $\frac{\epsilon \sum_{i=1}^N e^{y_i}}{\tau}$ replacing $e^y$. In other words, we just have to replace the constant volatility in the Black-Scholes equation by the average volatility during the time period under consideration.

This is precisely the substance of Merton's theorem. While we have assumed a specific process for the volatility, the astute reader will notice that the final result does not depend on that process as long as that process is independent of $x$.

An Extension to Merton's Theorem

We now extend Merton's theorem to the case of non-zero correlation for the stochastic process of the volatility that we are investigating. Since the present value of the option is given by

$$f(x, y, T) = \int_{-\infty}^\infty \matel{x, y}{e^{-\hat{H}\tau}}{x'} f(x', T) dx'$$ (109)

and $\matel{x, y}{e^{-\hat{H}\tau}}{x'}$ is given by equation (4.57) with $S_0$ and $S_1$ given by (4.59) and (4.51) respectively. Now, since $\int DY e^{S_0}$ describes the probability of a specific path for $y$ (we show this in detail in chapter 6), we see that the propagator can now be written as in equation (4.56) with $\omega,\, \theta,\, \eta$ and $\zeta$ being functionals of this path and $\nu$ being the final value of $V^{3/2-\alpha}$ for the path. Hence, if we generate paths for $y$ according to (4.6), the option price is given by

$$f(x, y, t) = \left\langle \int_{\ln K}^\infty dx' \frac{e^{S... ...u)}}{\sqrt{2\pi \epsilon (1-\rho^2)}} (e^{x'}-K) \right\rangle$$ (110)

(since the payoff of the option is given by $e^{x'}-K$) with the average taken over the paths for $V$. Since the propagator is in the form of a Gaussian, we can perform this integration to obtain

$$\begin{split} f(x, y, t) = &\left\langle SN(s_1)\exp\left( -... ...\zeta \right) - Ke^{-r\tau}N(s_2) \right\rangle \\ \end{split}$$ (111)

($N(x)$ denotes the cumulative normal distribution as in the previous chapter) where $V_i$ and $V_f$ are the initial and final volatilities of the path respectively and $s_1$ and $s_2$ are given by
$$\begin{align} s_1 = &\frac{\ln\left(\frac{S}{K}\right) + r\tau + \frac{1}{2} (1-... ...ta}{\sqrt{(1-\rho^2)\omega}}\\ s_2 = &s_1 - \sqrt{(1-\rho^2)\omega} \end{align}$$
The reader should be easily able to verify that equation (4.65) is the same as the Black-Scholes equation for any single volatility path when $\rho = 0$.

When $\alpha = \frac{1}{2}$, several simplifications occur. We see that $\theta = \zeta = \tau$ are known and $\eta = \omega$. In that case, (4.65) reduces to

$$\left\langle SN(s_1)\exp\left( -\left(\frac{\rho^2}{2} + \f... ...2}{2\xi}\right) \tau \right) - Ke^{-r\tau}N(s_2) \right\rangle$$ (112)

where $s_1$ and $s_2$ now have the relatively simple forms
$$\begin{align} s_1 = &\frac{\ln\left(\frac{S}{K}\right) + \left(r + \rho\left( \f... ...t)}{\sqrt{(1-\rho^2)\omega}}\\ s_2 = &s_1 - \sqrt{(1-\rho^2)\omega} \end{align}$$
Hence, we see that we have a straightforward solution for $\alpha = \frac{1}{2}$ even when the correlation is not zero.

When $\alpha = \frac{3}{2}$ and $\lambda = 0$, we obtain a similar simplification since $\zeta = \omega$ and $\eta = \tau$. In this case, we obtain the following expression for the option price

$$f(x, y, t) = \left\langle SN(s_1)\exp\left(-\frac{\rho}{2} \... ...\mu\tau}{\xi}\right) \right) - Ke^{-r\tau}N(s_2) \right\rangle$$ (113)

and $s_1$ and $s_2$ are now given by
$$\begin{align} s_1 = &\frac{\ln\left(\frac{S}{K}\right) + \left(r + \frac{\rho\mu... ...t)}{\sqrt{(1-\rho^2)\omega}}\\ s_2 = &s_1 - \sqrt{(1-\rho^2)\omega} \end{align}$$

For the case considered in Baaquie[8], we have $\lambda = 0$ and $\alpha = 1$. In this case, we have $\eta = \zeta = \int_0^\tau e^{y/2} dt$ which gives us

$$f(x, y, t) = \left\langle SN(s_1)\exp\left( -\frac{\rho^2\om... ...i^2}{2} \right) \eta \right) - Ke^{-r\tau}N(s_2) \right\rangle$$ (114)

where $s_1$ and $s_2$ are now given by
$$\begin{align} s_1 = &\frac{\ln\left(\frac{S}{K}\right) + r\tau + \left( -\rho^2 ... ...}{\sqrt{(1-\rho^2) \omega}}\\ s_2 = &s_1 - \sqrt{(1-\rho^2) \omega} \end{align}$$
which is somewhat more complicated since two functionals, $\omega$ and $\eta$ of the volatility path are involved. In this case, however, a perturbation analysis can be used to derive an approximate form for the probability distribution functions of the functionals. Due to this fortunate occurence, a series solution to this problem can be obtained.

The probability density function for the functionals is a very difficult quantity to obtain. The probability density function for $\int_0^\tau Vdt$ obtained for the special case $\alpha = 1/2$ in Stein and Stein[21]. Stein and Stein[21] have used this probability distribution function and the ``straightforward'' solution for $\rho = 0$ to get an analytic form of the solution for this case. We now see that the result can be extended to non-zero $\rho$ if we can find the joint probability density function of this functional and $V_i-V_f$. While the individual probability distribution functions can be obtained (the pdf for $\omega$ is obtained in Stein[21] and the pdf for $V_i-V_f$ is trivial), they are not independent.

Risk-Neutrality

We show that the expected value of the underlying security $S$ whose initial value is $S_0$ is given by $S_0e^{rt}$ after time $t$ has elapsed for a large class of stochastic processes including the one we are considering in this thesis. In other words, we show that $A = e^{-rt}S$ is a martingale.

We first change variables from $S$ to $A$ in (4.5) (changing $\phi$ to $r$ in accordance with risk-neutral valuation) to obtain

$$dA = Ae^{y/2} dz_1$$ (115)

(where $z_1$ is the time integral of $W$ and hence a Wiener process) where $y$ may depend on $A$ ($y$ can be stochastic). We now consider the more general process

$$dA = f(A, y) dz_1$$ (116)

We note that $E[dA] = 0$ so that $E[A(t+dt)\vert A(t) = A_0] = A_0$. In other words, we see that $A$ is a martingale. Hence, we have shown the result. In general, a martingale process cannot have a drift term.

While the result is simple, it has important consequences. We note that risk-neutrality alone cannot determine any constraints for the volatility process. Any volatility process whatsoever satisfies risk-neutrality.

Brownian Motion on a Riemmanian Manifold

We can also consider (4.5) and (4.6) as Brownian motion on a Riemmanian manifold. This way of treating stochastic differential equations is considered in detail in Zinn-Justin[30]. The metric tensor of the manifold for which (4.5) and (4.6) can be considered as simple Brownian motion is found to be

$$g = \frac{1}{\xi^2(1-\rho^2)} \left( \begin{matrix} \xi^2 e^... ...i e^{y(1/2-\alpha)} & e^{2y(1-\alpha)}\\ \end{matrix} \right)$$ (117)

Using this interpretation, we can derive the discrete form of the Lagrangian (4.37). Further, we can also derive a continuous time limit for the Lagrangian which is given by

$$\begin{split} L = \frac{1}{2} \left( (\dot{x}_i + f_i)g_{ij}... ...mbda e^{-x_2} - (1-\alpha)(\dot{x}_2 + f_2) \right) \end{split}$$ (118)

where $x_1 = x$, $x_2 = y$, $f_1 = r - \frac{e^y}{2}$, $f_2 = \lambda e^{-y} + \mu - \frac{\xi^2}{2}e^{2y(\alpha-1)}$ and the Einstein summation convention is implied. The reader will note that the continuous Lagrangian is the same as the limit of the discrete Lagrangian as $\epsilon \tendsto 0$ for the quadratic terms while the continuous Lagrangian contains some extra linear terms in $y$ and $\dot{y}$ which are related to the curvature of the manifold.

However, this Lagrangian is very difficult to deal with directly and hence, we use the discrete action to derive an efficient algorithm.