On Accurate and Provably Efficient
GARCH Option Pricing Algorithms
Yuh-Dauh Lyuu
∗Chi-Ning Wu
†Abstract
The GARCH model has been very successful in capturing the serial corre-lation of asset return volatilities. As a result, applying the model to options pricing attracts a lot of attention. However, previous tree-based GARCH op-tion pricing algorithms suffer from exponential running time, a cut-off maturity, inaccuracy, or some combination thereof. Specifically, this paper proves that the popular trinomial-tree option pricing algorithms of Ritchken and Trevor (1999) and Cakici and Topyan (2000) explode exponentially when the number of partitions per day, n, exceeds a threshold determined by the GARCH pa-rameters. Furthermore, when explosion happens, the tree cannot grow beyond a certain maturity date, making it unable to price derivatives with a longer maturity. As a result, the algorithms must be limited to using small n, which may have accuracy problems. The paper presents an alternative trinomial-tree GARCH option pricing algorithm. This algorithm provably does not have the short-maturity problem. Furthermore, the tree-size growth is guaranteed to be quadratic if n is less than a threshold easily determined by the model parameters. This level of efficiency makes the proposed algorithm practical. The surprising finding for the first time places a tree-based GARCH option pricing algorithm in the same complexity class as binomial trees under the Black-Scholes model. Extensive numerical evaluation is conducted to confirm the analytical results and the numerical accuracy of the proposed algorithm. Of independent interest is a simple and efficient technique to calculate the transition probabilities of a multinomial tree using generating functions.
Keywords: GARCH, trinomial tree, path dependency, option pricing, generating function
∗Corresponding author. Professor, Department of Computer Science & Information
Engi-neering, National Taiwan University, No. 1, Sec. 4, Roosevelt Rd., Taipei, Taiwan. E-mail: [email protected]. Also Professor of Department of Finance, National Taiwan University. The author was supported in part by NSC grant 92-2213-E-002-016.
1
Introduction
Efficient numerical algorithms play a critical role in derivatives pricing because it is of-ten imperative to obtain prices fast, particularly when prices change quickly. In both theory and practice, computational efficiency is measured in terms of running time. Two types of algorithms are usually distinguished: polynomial-time algorithms and exponential-time algorithms (Papadimitriou (1995)). Because the exponential func-tion grows extremely fast, exponential-time algorithms highlight the tradeoff between accuracy and speed much earlier than polynomial-time algorithms. Exponential-time algorithms are therefore said to suffer from combinatorial explosion and should be avoided wherever possible.
In the numerical pricing of derivatives, the continuous diffusion process for the asset price is often discretized to yield a tree first. Derivatives are then priced on the tree by standard backward induction. The lognormal diffusion, for instance, gives rise to the well-known CRR binomial tree of Cox, Ross, and Rubinstein (1979). Two critical features of the CRR tree, as well as its many binomial and trinomial variants, stand out: It recombines and an N-period tree contains only O(N2) nodes; i.e., it
exhibits a quadratic growth (see Fig. 1). As a consequence, simple derivatives such as vanilla options, barrier options, and lookback options can be efficiently priced as surveyed in Lyuu (2002). However, a polynomial-sized tree may still give rise to an exponential-time pricing algorithm if the derivative’s payoff is complex. The Asian option with a payoff depending on the arithmetic price average fits this characteri-zation. Specifically, the vast number of extra states (the running averages) needed by the Asian option’s path-dependent feature makes pricing on an N-period tree take time exponential in N. Approximations are mandatory for such derivatives to regain efficiency. Of course, approximation algorithms must not sacrifice accuracy. Important numerical techniques include the tree method with interpolation such as Hull and White (1993a), the PDE method such as Forsyth, et al. (2002), and the linear-programming technique such as Dempster and Richards (2000).
A fundamentally more difficult problem emerges when the explosion arises from the model itself. If the model generates trees that do not recombine, pricing is ex-pensive even for the simplest of derivatives like vanilla options. For example, when the volatility is not a constant, such as the interest rate model of Cox, Ingersoll, and Ross (1985), a brute-force discretization leads to exploding binomial trees that do not recombine. The problem may be rectified by the standard technique of Nelson and Ramaswamy (1990) that transforms the diffusion process into one with a constant volatility. But this methodology does not guarantee to do away with combinatorial explosion when the diffusion process is bivariate. For example, Chien (2003) demon-strated that the bivariate interest rate tree of Ritchken and Sankarasubramanian (1995) can explode exponentially even after the said transformation. The focus of
Figure 1: Binomial tree. Each node has two successor nodes. The number of nodes at any time t is t + 1, a linear growth in time. The total number of nodes of an
N-period binomial tree is PNt=0(t + 1) = (N + 2)(N + 1)/2 = O(N2), a quadratic
growth in maturity N.
this paper, the tremendously influential generalized autoregressive conditional het-eroskedastic (GARCH) model, is also bivariate.
Bollerslev (1986) and Taylor (1986) independently proposed the GARCH process for modeling the stochastic volatility of asset returns. Since then, the model has been generalized and used extensively in the finance literature on the modeling of financial time series; see Bollerslev et al. (1992) and Engle and Patton (2001). As the model has received strong empirical support, its application to option pricing draws a lot of attention. Duan (1995) was the first to propose a GARCH option pricing model. The massive path dependency of the pricing model initially favors Monte Carlo simulation over trees. But the Monte Carlo estimate is probabilistic; furthermore, options that can be exercised early, the so-called American options, can be accurately priced only with simulation schemes that employ advanced techniques. The appearance of the trinomial tree of Ritchken and Trevor (1999) addresses these problems and makes a strong case for trees. Their algorithm is simple and claims to be accurate and efficient. It is also general enough to work beyond GARCH models. GARCH option pricing techniques that are not based on trees include the Markov chain approximation of Duan and Simonato (2001), the Edgeworth tree approximation of Duan et al. (2002), and analytical approximations as in Heston and Nandi (2000). Among them, only the Markov chain approximation approach is capable of handling American options. This paper first investigates the performance of the Ritchken-Trevor algorithm and its modified version by Cakici and Topyan (2000). The findings are mixed, both theoretically and numerically. We prove that the Ritchken-Trevor-Cakici-Topyan (RTCT) algorithm can easily result in exponential-sized trees. This theoretical result is decisively backed up by numerical data. Interestingly, the simple condition for the combinatorial explosion mirrors that for GARCH to be nonstationary and is satisfied when the number of partitions per day, n, exceeds a typically small number. RTCT is hence efficient only if it is restricted to small n. But a small n does not guarantee
efficiency; in fact, even the smallest possible choice n = 1 can result in explosion under some circumstances. A small n can also yield inaccurate option prices. We note that the Cakici-Topyan algorithm is incompletely specified (see Section 4).
Even though raising n to improve accuracy will quickly incur exponential slow-down, one may be willing to trade efficiency for accuracy. (For example, the slow exponential-time algorithm can be run overnight on spare computers to obtain prices for situations likely to happen on the next trading day.) But this sensible practice turns out to be impossible for RTCT. Specifically, we prove that when explosion occurs, the RTCT tree cannot grow beyond a certain maturity, making it useless for pricing derivatives with a longer maturity. Numerical data again back up this analyt-ical result. The strong impossibility result obliterates the tradeoff between accuracy and efficiency for exploding RTCT trees. The loss of this tradeoff can be devastating if the permissible n’s fail to give accurate prices. The result also throws into question some of the calculated prices in Ritchken and Trevor (1999).
This paper then proposes a new trinomial-tree GARCH option pricing algorithm that addresses the above-mentioned shortcomings of RTCT. In a trinomial tree, the three successor nodes of each node must be such that their stock prices match the mean and variance of the model’s stock price asymptotically. The three branches must also carry valid branching probabilities. It is obvious that the higher the volatility, the wider the nodes must be spaced to meet the requirements. The RTCT tree, like typical trinomial trees, takes a flat middle branch from each node as in Fig. 2(a). Our new tree departs from that by making the middle branch track the expected stock price as closely as possible as in Fig. 2(b). Therefore, we call it the mean-tracking (MT) algorithm. By tracking the mean, the two flanking branches are expected to fan out less in their attempt to match the mean and variance of the next stock price. This in turn yields more compact trees. Although the mean-tracking idea works for any model, the effects are most striking for models with a stochastic volatility such as GARCH. The above argument is both intuitively appealing and provably true as this paper shows.
MT solves the short-maturity problem of RTCT. As a result, it accepts any
n without having to worry about the tree being cut short. The tradeoff between
efficiency and accuracy is thus restored. Because an exponential-sized tree consumes so many resources, one may want to be assured that they are not building one before committing the efforts. But is there a bound on n which guarantees a polynomial-sized MT tree? Indeed, such a bound exists and is surprisingly simple. Perhaps even more unexpectedly, the MT tree’s size is only quadratic in maturity if n does not exceed that bound. MT is therefore the first tree-based GARCH option pricing algorithm that is provably efficient. Interestingly, this result places MT in the same complexity class as CRR under the much simpler Black-Scholes model. A level of
(a) (b) (c)
Figure 2: Trinomial tree. (a) Three successor nodes follow each node. How widely the two flanking branches fan out around the middle branch depends critically on volatility. (b) The middle branch may maintain a drift to minimize that width. (c) When the vertical node spacing is a constant, the number of nodes at any time t is 2t + 1, a linear growth. The total number of nodes of an N-period trinomial tree is thus PNt=0(2t + 1) = (N + 1)2 = O(N2), a quadratic growth in maturity N.
efficiency like this opens up opportunities for the practical use of MT in pricing. The proof of the theoretical results requires only elementary techniques. Numerical experiments demonstrate that a small n gives accurate results. They also definitively confirm the quadratic node count.
In a general sense, the concept of mean tracking is not entirely new to the liter-ature. Numerical pricing techniques that remove the growth trend in the asset price can be said to adopt mean tracking, if only vaguely or implicitly at that. The idea, for instance, is explicit in the algorithms of Hull and White (1993b) and Li et al. (1995). Both works deal with the calibration of no-arbitrage interest rate models. The advantages of mean tracking in reducing complexity, however, have not been analytically explored until recently. For instance, Dai and Lyuu (2004) applied the mean-tracking idea to develop the first exact trinomitree Asian option pricing al-gorithm that breaks the long-standing 3N time barrier with a provable running time
of 3O(√N ), where N is the number of periods of the tree.
Because the numerical performance of the Cakici-Topyan (CT) version of the GARCH option pricing algorithm is slightly superior to that of the Ritchken-Trevor (RT) version, we will compare MT with CT in the numerical valuation of option prices. Our theoretical results on explosion and the shortened maturity of exploding trees apply to both RT and CT. Therefore, we will simply refer to both as RTCT in the analysis.
The paper is organized as follows. The GARCH model is presented in Section 2. Section 3 reviews the RTCT tree from which the MT tree derives. Several differences between the two will be pointed out along the way. A very simple and efficient way to calculate the transitional probabilities using generating functions is also presented there. Section 4 covers backward induction with interpolated volatilities, which are
needed to reduce the state space. In Section 5, a condition for the RTCT tree to explode exponentially is stated and proved. The nonexistence of the tradeoff between efficiency and accuracy for the exploding RTCT tree is proved in Section 6. Section 7 provides numerical data to back up the theoretical results and other negative claims regarding RTCT. Section 8 presents and analyzes MT. Section 9 proves that the MT tree grows only quadratically when n does not exceed a simple bound. Section 10 evaluates the various GARCH option pricing algorithms numerically. Section 11 concludes.
2
The GARCH model
Let St denote the asset price at date t and ht the conditional volatility of the
return over the (t + 1)th day [ t, t + 1 ]. Here, “one day” is just a convenient term for any elapsed time ∆t. The following risk-neutral process for the logarithmic price
yt ≡ ln St is due to Duan (1995): yt+1 = yt+ r − h2 t 2 + ht²t+1, (1) where h2t+1 = β0+ β1h2t + β2h2t(²t+1− c)2, (2)
²t+1 ∼ N(0, 1) given information at date t,
r = daily riskless return.
The model is bivariate as its state is described by (yt, h2t). Updating rule (2) for the
squared volatility, due to Engle and Ng (1993), is also called the nonlinear asymmetric GARCH or NGARCH for short. Other GARCH models are surveyed in Duan (1997). A positive c represents a negative correlation between the shock for the asset return and its conditional volatility. We assume β0, β1, β2 ≥ 0 to make the squared
volatilities h2
t positive. We further impose β1+β2 < 1 to make the model stationary.
The violation of a version of this inequality will be shown to make the RTCT tree explode. Specifically, explosion occurs when β1+β2n > 1 given c = r = 0 (the general
case of nonzero c and/or c is similar, but the explicit formula is complex). In contrast, the MT tree grows only quadratically in maturity when n ≤ (
q
1−β1
β2 − c)
2. Note that
when c = 0, this bound becomes β1 + β2n ≤ 1, the exact complement of the bound
for RTCT to explode. Hence both bounds are optimal in a sense. Throughout the paper, N will denote the maturity of the tree in days, which is also the maturity of the option to be priced by the tree.
3
The RTCT Tree
The RTCT trinomial tree approximates the continuous-state GARCH process with discrete states as follows. Each state is represented as a node. Partition a day into
n periods. Three states follow each state (yt, h2t) after a period. As the trinomial
tree recombines, 2n + 1 states at date t + 1 follow each state at date t. Let γ = h0
and γn = γ/
√
n . (Our theoretical results will turn out to be independent of how γ
is picked. Later, MT will choose a different γ for numerical reasons.) The tree is laid over nodes that are spaced by γn in their logarithmic prices as depicted in Fig. 3(a).
Consequently, the logarithmic price yt on each node equals y0+ kγn for some integer
k. (yt, h2t) 6 ?γn 6 ? ηγn -¾ 1 period -¾ 1 period (a) (b)
Figure 3: Jump parameter η and jump size ηγn. (a) The tree is laid over
a lattice whose vertically adjacent nodes are spaced by γn. (b) The two flanking
branches fan out around the middle branch to reach the two nodes that are η nodes away from the center. Although 2(η − 1) hollow nodes are not reached from the node on the left, they may be reached from other nodes. Here η = 3.
We next pick the jump size for state (yt, h2t). The jump size determines how much
the state’s three successor states are spaced. As emphasized earlier, the magnitude of the jump size depends on the volatility ht. By the geometry of the tree, the jump
size must be some integer multiple η of γn. We call η the jump parameter. The
jump parameter measures how much the two flanking branches fan out around the middle branch. It must be large enough for the three branches to match the mean and variance of yt+1. The three nodes one period hence extend over 2η + 1 nodes
(inclusively), and inductively the 2n + 1 nodes from (yt, h2t) at date t + 1 extend
over 2nη + 1 nodes (inclusively). See Fig. 3(b) for illustration. The larger the jump parameter η, the larger the tree because it extends over more nodes. The middle branch of the RTCT tree leaves the underlying asset’s price unchanged. In contrast, the MT tree lets the middle branch track the mean of yt+1, i.e., yt + r − (h2t/2).
This idea will result in a smaller jump parameter η, thus yielding a more compact tree. Figure 4 illustrates a 1-day trinomial tree with each day partitioned into n = 3 periods. (yt, h2t) 6 ? ηγn -¾1 period -¾ 1 day
Figure 4: RTCT trinomial tree for daily logarithmic price ytfor the duration
of one day. A day is partitioned into n = 3 periods, and the jump size is ηγn. The
7 values on the right approximate the distribution of yt+1 given (yt, h2t). Recall from
Fig. 3 that there are η − 1 nodes (which are not drawn), between any two vertically adjacent (black) nodes.
The probabilities for the up, middle, and down branches are
pu = h2 t 2η2γ2 + r − (h2 t/2) 2ηγ√n , (3) pm = 1 − h2 t η2γ2, (4) pd = h 2 t 2η2γ2 − r − (h2 t/2) 2ηγ√n . (5)
As the branching probabilities are picked to match the mean and variance of yt+1
given (yt, h2t) asymptotically, the tree converges to the continuous-state model (1).
From Eqs. (3)–(5), it is not hard to verify that valid branching probabilities exist (i.e., 0 ≤ pu, pm, pd≤ 1) if and only if
| r − (h2 t/2) | 2ηγ√n ≤ h2 t 2η2γ2 ≤ min µ 1 −| r − (h 2 t/2) | 2ηγ√n , 1 2 ¶ . (6)
The node count can be reduced by a factor of n by getting rid of the intraday nodes. The result is a (2n + 1)-nomial tree as in Fig. 5; it is multinomial with 2n + 1 branches from every state (yt, h2t). These 2n + 1 successor nodes extend over
2nη + 1 nodes (inclusively). This multinomial tree is the building block of the RTCT tree. Numerical data will demonstrate that although only 2n + 1 of those 2nη + 1 nodes are immediately reachable from (yt, h2t), the overwhelming majority of those
2nη + 1 nodes are reachable from the root state (y0, h20). The implication is that the
overwhelming majority of the nodes extended over by the RTCT tree are occupied by reachable states. (yt, h2t) A 6 ? ηγn -¾ 1 day
Figure 5: RTCT multinomial tree for daily logarithmic price yt for the
duration of one day. This heptanomial tree is the outcome of the trinomial tree in Fig. 4 after removing its intraday nodes. Recall that n = 3. In general, we infer from Fig. 3 that there are 2nη + 1 nodes at date t + 1 between the top and bottom nodes (inclusive) spaced γn apart, but only 2n + 1 of which are reachable from (yt, h2t) and
drawn above. The overwhelming majority of those 2n(η − 1) nodes not drawn are reached from the root state as will be demonstrated by numerical data later.
Updating rule (2) must be modified to reflect the adoption of the discrete-state tree model. State (yt, h2t) at date t is now followed by state (yt+ `ηγn, h2t+1) at date
t + 1, where h2 t+1 = β0+ β1h2t + β2h2t(²0t+1− c)2, (7) ²0 t+1 = `ηγn− (r − h2t/2) ht , ` = 0, ±1, ±2, . . . , ±n.
We will call the resulting tree the full RTCT tree. For example, node A in Fig. 5 contains state (yt+1, h2t+1), where
yt+1 = yt+ (−1)ηγn, h2 t+1 = β0+ β1h2t + β2h2t ½ (−1)ηγn− [ r − (h2t/2) ] ht − c ¾2 .
As part of a larger RTCT tree, node A may contain other states which are reached from states other than (yt, h2t).
From the underlying trinomial model, the transition from state (yt, h2t) to state
(yt+ `ηγn, h2t+1) happens with probability
P (`) ≡ X ju,jm,jd n! ju! jm! jd! pju u pjmmpjdd, (8)
where ju, jm, jd ≥ 0, n = ju + jm + jd, and ` = ju − jd. This seemingly
compli-cated formula for probabilities P (`) can be calculated very efficiently with the simple generating-function technique in Lyuu (2002) as follows. Note that
(pux + pm+ pdx−1)n= n
X
`=−n
P (`) x`.
Therefore, we can multiply out (pux+pm+pdx−1)nand then retrieve the probabilities
P (`) by reading off the coefficients of the powers of x. The computation takes O(n2)
steps. Compared with the complex formula (8), the generating-function approach is straightforward and stable.
As volatility ht changes through time, we may have to pick different jump
pa-rameters η for different states so that all branching probabilities pu, pm, and pd
lie between 0 and 1. This entails varying jump sizes. As the necessary requirement
pm ≥ 0 implies
η ≥ ht/γ, (9)
RTCT goes through
η = d ht/γ e, d ht/γ e + 1, d ht/γ e + 2, . . . (10)
until valid probabilities are obtained or until their nonexistence is confirmed by in-equalities (6). The latter case means the tree cannot grow further. Later proofs depend critically on noting that the magnitude of η grows at least as fast as ht by
virtue of inequality (9). Hence when ht grows exponentially, the resulting tree must
do likewise.
Every node at date t on the tree carries a different logarithmic price yt. However,
more than one path from the root state (y0, h20) may lead to the same node, each
yielding a different squared volatility h2
t. The number of possible values of h2t at a
node thus equals the number of paths reaching the node. Each h2
t picks its own jump
parameter η. Figure 6 depicts a 3-day tree with n = 1. Nodes A and B each have one
h2
t, whereas node C has two. The h2t at nodes A and B pick the same jump parameter
η = 2, whereas the two h2
t at node C pick different jump parameters η = 1, 2. The
overall tree structure is irregular because of the varying jump parameters. Hollow nodes in Fig. 6 are not occupied because they are unreachable from the root state
(y0, h20). As will be shown later, their count is minuscule. This outcome is expected
by observing node C. Although the three branches from node A do not reach node C, node C is nevertheless reached from the two nodes below node A. Another factor contributing to making unreachable nodes few is the fact that most nodes contain more than one squared volatility h2
t. If these squared volatilities pick different jump
parameters such as those at nodes C and D, more than three branches will emanate from the node. More branches make more nodes reachable.
(y0, h20) A B C D 6 ? γn = γ1 -¾ 3 days
Figure 6: Geometry of a 3-day RTCT tree. A day is partitioned into n = 1 period. Nodes A and B have a jump parameter of 2. Nodes C with two h2
t and D
with three h2
t have two jump parameters: 1 and 2. All other nodes have a jump
parameter of 1. Hollow nodes are not reachable from the root state (y0, h20).
4
Interpolated Volatilities and Backward
Induc-tion
The number of possible volatilities ht at a node equals the number of paths reaching
it. An algorithm that keeps all of these volatilities cannot be efficient as their count is exponential. Therefore, the full RTCT tree (7) must be approximated. The standard approximation methodology by Hull and White (1993a) and Ritchken, Sankarasubra-manian, and Vijh (1993) is adopted by both RTCT and, later, MT. The actual RTCT tree now keeps only the maximum volatility hmax and the minimum volatility hmin
hmax .. . ht .. . hmin .. . ht+1 .. .
-Figure 7: Case where maximum volatility follows an interpolated volatility. The maximum volatility ht+1 at the node on the right follows interpolated volatility
ht.
at each node. It then creates K − 2 volatilities between hmin and hmax (what values
these volatilities assume will be elaborated soon). Let us call these K − 2 volatilities interpolated volatilities because they are not the results of actually following updat-ing rule (7). Instead, they are artificial volatilities generated via interpolation. The number of volatilities per node is now a constant K.
CT and RT differ in the way they calculate hmin and hmax at each node. RT
calculates them with the interpolated volatilities taken into consideration. This means that every node at date i generates 2n + 1 volatilities for each of its K volatilities via the updating rule. The number of volatilities generated per node is hence K(2n + 1). These volatilities together with their associated branches then determine the hmin and hmax at every node at date i + 1. CT, on the other hand, calculates hmin and hmax in
the tree-building process without the interpolated volatilities. This means that every node at date i generates 2n+1 volatilities for each of its two volatilities hminand hmax
via the updating rule. The number of volatilities generated per node is hence only 2(2n + 1). These volatilities together with their associated branches then determine the hmin and hmax at every node at date i + 1. Both CT and RT use interpolated
volatilities in backward induction.
The interpolated volatilities between hmin and hmaxat a node are artificial because
they are not the results of applying updating rule (7) starting from the root state (y0, h20) at time 0. For the same reason, hmin and hmax in CT are not artificial as they
are the results of following chains of updating rule without interpolation. For RT, even hmin and hmax may be artificial if they are the results of applying the updating
rule to interpolated volatilities of the previous date (see Fig. 7 for illustration). In RTCT, the K squared volatilities at a node are equally spaced between h2
min and h2 max: h2 min+ j h2 max− h2min K − 1 , j = 0, 1, 2, . . . , K − 1.
hmax .. . ht .. . hmin ¾ .. . ht+1 .. .
Figure 8: Backward induction. Volatility ht+1 follows ht by the updating rule.
Because ht+1 does not match any interpolated volatility, its corresponding option
value is found by interpolating from the two option values whose volatilities bracket it.
To be specific, the K logarithms of squared volatilities are equally spaced between ln h2
min and ln h2max in MT; they are
exp · ln h2min+ j ln h 2 max− ln h2min K − 1 ¸ , j = 0, 1, 2, . . . , K − 1.
Smaller volatilities are thus sampled more finely than larger volatilities. We call it the log-linear interpolation scheme.
After the tree is built, backward induction commences. For the volatility ht+1
following state (yt, h2t) via updating rule (7), RTCT finds the two volatilities that
bracket ht+1. The option price corresponding to ht+1 is then interpolated linearly
from the option prices corresponding to the bracketing volatilities. Figure 8 illustrates this procedure for a branch. After the option prices from all 2n + 1 branches are available, the option price for state (yt, h2t) is calculated as their average discounted
value weighted by the branching probabilities.
Before closing this section, we remark that CT is under-specified. CT maintains only hmin and hmax at each node in growing the tree. The interpolated volatilities at
each node are not involved in the tree-building process. Only in backward induction are the K − 2 interpolated volatilities added to approximate the full RTCT tree. It is therefore possible for an interpolated volatility’s successor volatility to reach a node that is not reached during the tree-building process and thus has no option prices at all. Take Fig. 8 for example. If ht+1 lies in a node not reached during
tree building, the node will not have option values at all and backward induction must be terminated. Rare as these situations are, they do arise when n and N are both large. RT, in contrast, has no such problems because the same volatilities (interpolated or otherwise) used in building the tree are used in backward induction. In other words, the interpolated volatilities’ branches have all been followed through in the tree-building process. MT follows RT in using all interpolated volatilities, not just hmin and hmax, in growing the tree.
5
Sufficient Condition for RTCT To Explode
It is common practice to raise numerical accuracy by increasing n. Unfortunately, this section proves that the largest value of ht at date t grows exponentially in t if n
exceeds a threshold. When this happens, the value of η must also grow exponentially by virtue of inequality (9). Note that the 2n + 1 successor nodes of a state extend over 2nη + 1 nodes (recall Fig. 5). So when η grows exponentially, the RTCT tree explodes. We conclude that the RTCT tree must be restricted to small n to have any hope of being efficient. A small n is no guarantee that the tree will not explode, however, because our threshold for explosion is a sufficient condition only. Hence it does not follow that an n not exceeding that threshold escapes explosion. Because RTCT employs a search procedure (10), finding the necessary condition for explosion seems difficult.
We now provide the simple argument for the exponential growth of the largest value of ht at date t for RTCT. Assume r = 0 and c = 0. Updating rule (7) is now
h2 t+1 = β0+ β1h2t + β2 £ `ηγn+ (h2t/2) ¤2 , ` = 0, ±1, ±2, . . . , ±n. Set ` = n to make h2
t+1 as large as possible. The updating rule now becomes
h2 t+1 = β0+ β1h2t + β2 £ √ n ηγ + (h2 t/2) ¤2 as γn = γ/ √ n ≥ β0+ β1h2t + β2 £ √ n ht+ (h2t/2) ¤2 as ηγ ≥ ht ≥ β0+ β1h2t + β2nh2t = β0+ (β1+ β2n) h2t. By induction, h2t+1 ≥ β0 t X i=0 (β1+ β2n)i+ (β1 + β2n)t+1h20 = β0 1 − (β1+ β2n) + · h2 0+ β0 (β1+ β2n) − 1 ¸ (β1+ β2n)t+1.
The above expression grows exponentially in t if
β1+ β2n > 1. (11)
Interestingly, this inequality is reminiscent of the necessary condition β1+ β2 ≥ 1 for
GARCH to be nonstationary.
When r 6= 0 or c 6= 0, the largest value of ht at date t still grows exponentially
in t as long as n is suitably large. The argument is more tedious but essentially identical. We conclude that the RTCT tree grows exponentially if n is large enough.
6
The Shallowness of the Exploding RTCT Tree
Can a large n be picked assuming we are willing to accept long running times? Unfortunately, as will be shown below, RTCT does not admit such a tradeoff between accuracy and speed because it will be forced to be short-dated when explosion occurs. Consequently, even if one is willing to accept inefficiency, the RTCT tree may not grow to the desired maturity date. An exploding RTCT tree is therefore of very limited use in practice for reasons other than exponential complexity.
The technical reason for shallowness is that there is a global ceiling on volatility
ht for valid branching probabilities to exist: h2t ≤ 2(r + n) + 2
√
2rn + n2. This upper
bound depends only on r and n. So when the largest value of ht grows exponentially
in t, this ceiling will be quickly reached at a small t and the tree can grow no further. We next derive the promised upper bound that leads to the impossibility result.
Inequalities (6) imply | (h2 t/2) − r | 2ηγ√n ≤ h2 t 2η2γ2, h2 t 2η2γ2 ≤ 1 2. Hence h2 t ≤ (ηγ)2 ≤ · h2 t √ n | (h2 t/2) − r | ¸2 ,
which can be simplified to yield £
(h2t/2) − r¤2 ≤ nh2t.
The above quadratic inequality (in h2
t) is equivalent to 2(r + n) − 2√2rn + n2 ≤ h2 t ≤ 2(r + n) + 2 √ 2rn + n2. We conclude that h2 t ≤ 2(r + n) + 2 √ 2rn + n2 (12)
is necessary for the existence of valid branching probabilities for RTCT. As promised earlier, the bound on h2
t does not depend on the choice of γ as the equation γ = h0
never enters the analysis.
This impossibility finding may sound puzzling at first. After all, under the Black-Scholes model, valid branching probabilities always exist as long as n is suitably large. Why, one may ask, can’t the same property hold here? The answer lies in the volatility process. The daily volatility in the Black-Scholes model is a constant, which amounts to setting ht to some fixed number. So every state solves the same
Eqs. (3)–(5) for the probabilities, and increasing n will eventually have inequality (12) satisfied for all states. In contrast, the volatility fluctuates under GARCH. Thus each
state (yt, h2t) faces different Eqs. (3)–(5) in solving for the probabilities. Increasing
n beyond a certain threshold makes it impossible to satisfy inequality (12) when t is
suitably large. The reason is that the largest h2
t grows exponentially in t, eventually
breaching the upper bound at a suitably large t.
7
Numerical Evaluation of RTCT
The following parameters from Ritchken and Trevor (1999) and Cakici and Topyan (2000) will be assumed throughout the section unless stated otherwise: S0 = 100, r = 0, h2
0 = 0.0001096, γ = 0.010469, β0 = 0.000006575, β1 = 0.9, β2 = 0.04, and c = 0. As r = c = 0, criterion (11) says combinatorial explosion occurs when
n > 1 − 0.9
0.04 = 2.5.
Figure 9 confirms the exponential growth of the RTCT tree with n = 3, 4, 5 . Note that the total node count increases with n. For comparison, the standard trinomial tree contains only 2t + 1 nodes at date t.
25
50
75 100 125 150 175
Date
5000
10000
15000
20000
25000
Figure 9: Exponential growth of the RTCT tree. Each curve plots the number of nodes at date t of the tree. The dotted line is based on n = 3, the dashed line on n = 4, and the solid line on n = 5. The parameters are S0 = 100, r = 0, h2
0 = 0.0001096, γ = 0.010469, β0 = 0.000006575, β1 = 0.9, β2 = 0.04, and c = 0.
The number of nodes (call it M) is critical because the running time is proportional to KM . We mention earlier that there may be nodes extended over by the tree which are not reachable and are depicted as hollow nodes in Fig. 6. In principle, if such nodes are numerous, RTCT can potentially run more efficiently by adopting clever
programming techniques to skip them. But Figure 10 shows that the proportion of unreachable nodes is small for n = 3, 4, 5. We will see shortly that the same conclusion holds for larger n as well. As the overwhelming majority of nodes are reachable, skipping unreachable ones brings no substantial benefits. Therefore, by tree size we shall mean the number of nodes extended over by the tree regardless of reachability.
25
50
75
100 125 150 175
Date
0.5
1
1.5
Figure 10: The percent of unreachable nodes. The plots show the percent of unreachable nodes among all nodes at each date. The dotted line is based on
n = 3, the dashed line on n = 4, and the solid line on n = 5. The number of
unreachable nodes is insignificant in all 3 lines. The parameters are S0 = 100, r = 0, h2
0 = 0.0001096, γ = 0.010469, β0 = 0.000006575, β1 = 0.9, β2 = 0.04, and c = 0.
Suppose we select n = 100 to seek very high accuracy at the expense of efficiency. The theory predicts that the RTCT tree’s final maturity will be cut short. Indeed, with r = 0, inequality (12) imposes the universal upper bound of h2
t ≤ 4n = 400.
This means a node with ht> 20 cannot have valid branching probabilities and thus
the tree cannot grow beyond date t. As this ceiling is breached somewhere at date 9 because of the exponential growth of the largest value of ht, the tree stops growing
there. Table 1 lists the final dates under various n. Observe that the tree’s final maturity date decreases rapidly as n increases. For example, it is 72 with n = 5, 34 with n = 10, and 12 with n = 50. To be useful, n cannot be so large as to make the tree’s final maturity date fall short of the derivative’s. This concern sets a firm upper bound on n. Table 1 also tabulates the total numbers of nodes and unreachable ones among them. Again, the overwhelming majority of the nodes are reachable.
Some of the calculated option prices in Ritchken and Trevor (1999) use n as large as 25 and option maturity dates as far as 200. These choices contradict our analysis
and data. For example, according to Table 1, the RTCT tree with n = 25 stops growing as early as date 18. Therefore, some of their prices should be viewed with caution.
Cakici and Topyan (2000) used n = 1 throughout their paper. Explosion is avoided for their choice of parameters, which are also used in this section. But there are two problems with this choice of n. First, CT does not always generate accurate option prices with n = 1. For example, suppose we select β0 = 0.000007 instead of
0.000006575. The option prices in Table 2 show that CT deviates from the Monte Carlo estimates. Second, it is not true that the choice n = 1 always escapes explosion. To support this point, we use the same parameters as Fig. 9 except for c = 2 and
n = 1. The resulting CT tree explodes as plotted in Fig. 11.
10
20
30
40
50
Date
500
1000
1500
2000
2500
3000
Figure 11: Exponential growth of the CT tree with n = 1. This CT tree explodes. The final maturity date of the tree is 54. The parameters are S0 = 100, r = 0, h2
0 = 0.0001096, γ = 0.010469, β0 = 0.000006575, β1 = 0.9, β2 = 0.04, c = 2,
and n = 1.
One of the most critical questions to ask of an approximation algorithm such as RTCT is whether it converges with increasing n. Table 3 and Figure 12 show that there is a downward trend in the calculated option prices except for the very short maturity of N = 2 days. Moreover, the downward trend accelerates as n increases. To rule out the possibility that the problem originates from the discrete-state full RTCT tree approximation of the continuous-state model, Monte Carlo simulation of the full RTCT tree is carried out. The data in Table 3 demonstrate that the tree model produces Monte Carlo estimates generally consistent with those of the continuous-state model. Hence the problem lies mainly with the way the full RTCT tree is in turn approximated. We will see later that a major reason causing RTCT
to deviate numerically from the true option price is its use of the linear volatility interpolation scheme. ! "# $% & $''!" # $% & ( ) *+, -./ 0 1 2 3 4 56 7 8 9 :; !" # $% & $''!" # $% & <= > 0 0 0 0 0? 0 0 0> ? 0 @A B@C DEFGHDIJ K B@ I LL FGHDIJ K ( ) M +, -./ 0 0 0? 0 00 ? ? N O P Q RS T U V W X YZ [Y\ ]^_` a ]bc d [ Y bee_` a ]bc d
Figure 12: Select option prices from Table 3. MC lower bound equals ∞L; MC
upper bound equals ∞U.
8
The Mean-Tracking (MT) Tree
To recap, RTCT has at least four weaknesses. First, it explodes exponentially when
n exceeds a threshold. Second, it is not known whether there is a simple formula
for the threshold n∗ so that the RTCT tree escapes explosion as long as n ≤ n∗.
Third, when explosion happens, the tree is cut short, making it unable to price derivatives with a longer maturity date. Fourth, option prices may fail to converge as n increases. We next turn to MT that addresses these problems. MT makes two critical changes to RTCT. The first is to replace the linear interpolation scheme
with the log-linear interpolation scheme. This addresses the convergence problem mentioned earlier. The second is to let the middle branch of the multinomial tree track the mean of yt+1. This addresses the explosion problem and its consequence
of shortened maturity. It also yields a simple formula for the above-mentioned n∗.
Surprisingly, MT’s tree size grows only quadratically when n ≤ n∗. This result makes
the MT tree size asymptotically the same as that of the CRR binomial tree. The analysis of MT will turn out to require only elementary techniques.
8.1
Volatility Interpolation Schemes
The distribution of the volatilities reaching a node plays a key role in pricing ac-curacy. RTCT essentially assumes that the distribution is uniform: Interpolated squared volatilities are equally spaced between the minimum and the maximum squared volatilities. Figure 13, however, shows that the actual distribution is closer to a lognormal distribution than a uniform one. It strongly suggests that there be more interpolated volatilities at the lower values than at the higher values. This is the rationale for MT’s adopting the log-linear interpolation scheme, in which the logarith-mic volatilities are equally spaced. The log-linear scheme is also used by the Markov chain approximation of Duan and Simonato (2001) for the same numerical consid-erations. Similar findings in the case of the Asian option exist. For example, Dai (2004) proved that linear interpolation schemes result in overestimates, and Forsyth, et al. (2002) demonstrated that linear interpolation schemes may not converge to the correct price.
8.2
Tree Building
At date t, let node A be the node closest to the mean of yt+1 given (yt, h2t), which
equals yt+ r − (h2t/2). For convenience, we use µ to denote this conditional mean
minus the current logarithmic price:
µ ≡ r − h2t
2
(see Fig. 14). By the geometry of the tree, node A’s logarithmic price equals yt+ aγn
for some integer a. The criterion by which node A is chosen ensures that
| aγn− µ | ≤
γn
2 . (13)
To create the desired multinomial tree, make the middle branch of the (2n + 1)-nomial tree line up with node A as in Fig. 15. Although a node reaches only 2n + 1 nodes after one day, the top and bottom nodes extend over
Figure 13: Histogram of volatility distribution. All parameters are from Table 3 with n = 10 and N = 10. Monte Carlo simulation of the full RTCT tree (7) is used to record the volatilities at the middle node at maturity. Out of 500,000 paths, about 10,000 reach the node. The recorded minimum and maximum squared volatilities are 0.000085 and 0.000217, respectively. Squared volatilities h2
t are multiplied by
1,000,000 for easier reading.
nodes as in RTCT. Amazingly, as will be demonstrated later, the overall node count is only quadratic in N if n does not exceed a threshold. The probabilities for the upward, middle, and downward branches equal
pu = nh2 t + (aγn− µ)2 2n2η2γ2 n − aγn− µ 2nηγn , pm = 1 − nh2 t + (aγn− µ)2 n2η2γ2 n , pd = nh2 t + (aγn− µ)2 2n2η2γ2 n + aγn− µ 2nηγn .
As they match the mean and variance of the GARCH process at date t + 1 asymp-totically, convergence is guaranteed. State (yt, h2t) at date t is followed by state
γ γ 1 day
µ
Figure 14: The next middle node via mean tracking. The cross identifies the true mean of yt+1. Two nodes, A and B, bracket it. Between them, node A has a
logarithmic price closer to the mean. The number aγn denotes the difference between
yt and node A’s logarithmic price.
(yt+ `ηγn, h2t+1) at date t + 1, where h2t+1 = β0+ β1h2t + β2h2t(²00t+1− c)2, (15) ²00t+1 = `ηγn+ aγn− (r − h 2 t/2) ht , ` = 0, ±1, ±2, . . . , ±n.
We will call the resulting tree the full MT tree. Formulas for the transition probabil-ities are stated in equation (8). As the full MT tree contains an exponential number of states, it must be approximated. But unlike in the RTCT case, the approximated MT will be efficient, accurate, and never cut short.
The conditions for the probabilities to lie within 0 and 1, i.e., 0 ≤ pu, pm, pd≤ 1,
are | aγn− µ | 2nηγn ≤ nh 2 t + (aγn− µ)2 2n2η2γ2 n , (16) nh2 t + (aγn− µ)2 2n2η2γ2 n ≤ 1 2, (17) nh2 t + (aγn− µ)2 2n2η2γ2 n ≤ 1 − | aγn− µ | 2nηγn . (18)
n
a
γ
nηγ
)
,
(
y
th
t2 nηγ
Figure 15: The MT multinomial tree for daily logarithmic price yt for a
duration of one day. A day is partitioned into n = 3 periods, and the three jump sizes in each period are (aγn/n)+ηγn (upward), aγn/n (middle), and (aγn/n)−ηγn
(downward). The central branch of the tree lines up with node A, the node closest to the mean of yt+1 as stated in Fig. 14. The solid nodes are actually used in pricing,
but the gray nodes are for illustration only. This heptanomial MT tree should be compared with the RTCT counterpart in Fig. 5.
Inequalities (16)–(17) are equivalent to p nh2 t + (aγn− µ)2 nγn ≤ η ≤ nh 2 t + (aγn− µ)2 nγn| aγn− µ | . (19)
Inequalities (13) and (17) together imply inequality (18) because
nh2 t + (aγn− µ)2 2n2η2γ2 n + | aγn− µ | 2nηγn ≤ 1 2 + 1 4nη ≤ 1.
Hence the probabilities are valid if and only if the much simpler inequalities (19) hold. MT can avoid the short-maturity problem of RTCT. Let H2
min ≡ min(h20, β0/(1 − β1)) to make Hmin2 ≤ h2t for t ≥ 0 (see Appendix A for justification). Now, a valid
integer η can always be found regardless of the value of n as long as γn is less than
some constant. More precisely, interval (19) contains positive integers for the jump parameter η to take its value in when
γn2 ≤ Hmin2 (20) as proved in Appendix B. With the existence of η guaranteed, MT will never be cut short, as promised.
Rather than searching for an η to satisfy inequalities (19), MT simply sets η = &p nh2 t + (aγn− µ)2 nγn ' . (21)
Although other choices are clearly possible, this particular choice is amenable to the analysis on the size of the MT tree later.
We proceed to choose γ, hence γn as well because γn = γ/
√
n . If γ ≤ Hmin,
then γn satisfies inequality (20) for all n. A smaller γ generally leads to larger trees,
hence longer running times. On the other hand, a smaller γ is expected to result in better accuracy because of the finer grain. To strike an overall balance between accuracy and convergence speed, we set γ = Hmin/2; thus
γn=
Hmin
2√n. (22)
All the parameters of MT have now been specified.
9
A Sufficient Condition for the Nonexplosion of
MT
In practice, it is essential to know beforehand that the chosen n will not result in an exploding tree before tree building is even attempted. Without this knowledge, tree building may take a long time if the tree turns out to explode and may even end up with a tree not meeting the required maturity if shortened maturity is an issue. In the case of MT, the criterion for nonexplosion is simple: The MT tree does not explode if n ≤ Ãs 1 − β1 β2 − c !2 . (23)
Note that the threshold is independent of r. When inequality (23) holds, the tree-size growth is only quadratic in maturity, the same as the CRR tree. (See Appendix C for proof.) The MT tree is thus the first tree-based GARCH option pricing algorithm that is provably efficient. This surprising finding makes MT a very practical algorithm.
Earlier, we show that the RTCT tree explodes if n exceeds some threshold; in fact, the threshold is (1 − β1)/β2 when c = r = 0. But we are silent about the
tree size when that threshold is not breached. The positive result (23) for MT fills that void because it says the MT tree is efficient if n does not exceed an explicit threshold. Surprisingly, the sufficient condition for MT’s nonexplosion reduces to
n ≤ (1 − β1)/β2 when c = 0. As the same (1 − β1)/β2 is the threshold of explosion
10
Numerical Evaluation of MT
RTCT can be inaccurate in pricing options. This fact is documented in Table 4. In every combination of n and K, RTCT deviates from the simulation results. In contrast, MT produces prices within the 95% confidence interval for K ≥ 2 with
n = 1, K ≥ 10 with n = 2, and K ≥ 50 with n = 3, 4. Hence MT succeeds
where RTCT fails. Observe that K needs to be increased for a larger n because the resulting increase in the number of volatilities per node demands more resolution.
n
! " # $ % ! " #$ % ! " #$ % ! " # $ % #&&! " #$ % # && !" # $ % # && !" # $ % #&&! " #$ % ' ' ' ' ( )*+ ,-. ( )*+ ,-. ( )*+ ,-. ( )*+ ,-. / / / /n
0 1 2 3 45 6 7 8 9: 0 1 2 3 45 6 7 8 9: 0 1 2 3 45 6 7 8 9: 0 1 2 3 45 6 7 8 9: ! " #$ % ! " #$ % ! " #$ % ! " #$ % #&&! " #$ % #&&! " #$ % #&&! " #$ % #&&! " #$ % ; ;; ; < =>? @AB < =>? @ A B < =>? @ A B < =>? @AB C C C C / / / / / / / / / / / / / / / / /D /D /D /D / / / / / / / / / C /C /C / C DDDD ////n
0 1 2 3 45 6 7 8 9: 0 1 2 3 45 6 7 8 9: 0 1 2 3 45 6 7 8 9: 0 1 2 3 45 6 7 8 9: EF EF EF EF GEH IJKL M INO P GEH IJK LM INO P GEH IJK LM INO P GEH IJKL M INO P GE NQQKL M INO P GE NQQKL M INO P GE NQQKL M INO P GE NQQKL M INO P GF GF GF GF R R R R S TUV WXY S TUV WXY S TUV WXY S TUV WXY / / / / / / / / /D /D /D /D / / / / // // // // D D D D DDDDn
Z [ \ ] ^_ ` a b c d Z [ \ ] ^_ ` a b cd Z [ \ ] ^_ ` a b cd Z [ \ ] ^_ ` a b c d ef e f e f ef g eh ijklmino p geh ijklmino p geh ijklmino p g eh ijklmino p g e nqqklmino p ge nqqklmino p ge nqqklmino p g e nqqklmino p gf gf gf gfFigure 16: Select option prices from Table 5. MC lower bound equals ∞L; MC
upper bound equals ∞U. The option prices under the CT algorithm come from
Fig. 12.
Next we benchmark MT’s performance with increasing n. For this purpose, we duplicate the settings in Table 3 for RTCT to produce Table 5 for MT. Select prices are plotted in Fig. 16 for graphical depiction. Unlike RTCT, all prices generated by MT are within the 95% confidence interval of Monte Carlo simulation of the full MT
tree (15). Furthermore, MT provides results very close to the true option price even with a small n. MT therefore outperforms RTCT in terms of overall accuracy and convergence speed.
All the numerical experiments up to now assume r = c = 0. In Table 6, we investigate MT’s options pricing accuracy with nonzero r and c, specifically, r = 5% (annual) and c = 0.5 from Duan and Simonato (2001). Although a few of the computed option prices are outside the 95% confidence interval, they are nontheless quite close to the Monte Carlo estimates. Furthermore, they are as good as those prices in Table 3 of Duan and Simonato (2001) that are allowed the most computation times.
Criterion (23) is a sufficient condition for MT to be efficient. In this case, the number of tree nodes at date t is the linear O(t); thus the total node count becomes quadratic in maturity. We next use two concrete cases to cross-check the theoretical result. The first setting adopts β1 = 0.8, β2 = 0.1, and c = 0. The tree should not
explode because 1 = n < Ãr 1 − 0.8 0.1 − 0 !2 = 2.
Indeed, the tree size grows linearly with date as shown in Fig. 17. The total number of nodes is therefore quadratic in maturity, in complete agreement with the theory. Take another setting with β1 = 0.8, β2 = 0.1, and c = 0.9. Criterion (23) is now
violated because 1 = n > Ãr 1 − 0.8 0.1 − 0.9 !2 = 0.264416.
The tree also turns out to grow exponentially as shown in Fig. 18. Unlike the RTCT tree, however, an exploding MT tree remains useful because it will not be cut short. Criterion (23) is a sufficient, but not a necessary, condition for the MT tree to escape explosion. Hence its violation does not necessarily imply exponential running time. Take the parameters β1 = 0.8, β2 = 0.1, and c = 0.5 in Table 6 for example.
The criterion for nonexplosion is violated for all n ≥ 1. But in fact, the case of n = 1 does not result in a tree of exponential size even though the total node count is more than quadratic. Explosion sets in for n ≥ 2.
11
Conclusions
GARCH option pricing is difficult because of the GARCH model’s bivariate and path-dependent nature. We prove that the Ritchken-Trevor-Cakici-Topyan (RTCT) GARCH option pricing algorithm is inefficient as the RTCT tree explodes for n ex-ceeding a threshold. Worse, a global upper bound on the volatility renders the tree shallow when explosion occurs. Cakici and Topyan (2000) claimed that RTCT is
10 20 30 40 50 60 70 t 100
200 300 400
Figure 17: Linear growth of the MT tree. The parameters are from Table 6 except for n = 1, c = 0, N = 74, and K = 20. Because this set of parameters satisfy criterion (23), the tree should not explode. Indeed, the number of nodes at date t grows linearly with t; the total number of nodes hence grows quadratically in maturity N. 10 20 30 40 50 60 70 t 5000 10000 15000 20000
Figure 18: Exponential growth of the MT tree. The parameters are identical to those in Fig. 17 except for c = 0.9. This set of parameters violate criterion (23). The number of nodes grows exponentially with date t.
empirically accurate with n = 1 for vanilla options. Unfortunately, numerical data demonstrate that both inaccuracy and even explosion can result with such under-refined trees. These theoretical results literally carry over to the BDT-GARCH in-terest rate model of Bali (1999).
The paper modifies RTCT to obtain a simple tree, the mean-tracking (MT) tree. The MT tree is both accurate and provably efficient when n does not exceed a simple threshold. Surprisingly, the tree-size growth is only quadratic in maturity if n does not exceed the threshold. This is the first tree-based GARCH option pricing algorithm that provably does not explode if certain conditions are met. The complexity in fact matches that of the influential and popular CRR binomial tree. Both the threshold and the quadratic size are optimal. Finally, the MT tree does not
suffer from the short-maturity problem of the RTCT tree. We conclude that MT is a provably efficient algorithm for derivatives pricing under the GARCH option pricing model. All our theoretical results are proved and backed up by extensive numerical experiments.
Appendix A
Updating rule (15) says h2
t ≥ β0 + β1h2t−1. By applying the inequality repeatedly,
h2t ≥ β0+ β1h2t−1 ≥ · · · ≥ β0(1 + β1+ β12+ · · · + β1t−1) + β1th20 = β0 1 − β1 + βt 1 µ h2 0− β0 1 − β1 ¶ . (24) Suppose h2
0 ≥ β0/(1 − β1). Then inequality (24) implies h2t ≥ β0/(1 − β1) for t ≥ 0.
On the other hand, suppose h2
0 < β0/(1 − β1). Let h20 = β0/(1 − β1) − ² for some ² > 0. By inequality (24) again, h2 t ≥ β0 1 − β1 − βt 1² ≥ β0 1 − β1 − ² = h2 0
because β1 < 1. So h2t ≥ h20 for t ≥ 0. Combine the two cases to yield h2
t ≥ min(h20, β0/(1 − β1)), t ≥ 0.
Appendix B
To confirm that interval (19) contains an integer, it suffices to establish p nh2 t + (aγn− µ)2 nγn ≤ nh2t + (aγn− µ)2 nγn| aγn− µ | − 1. (25) First, nh2 t + (aγn− µ)2 nγn| aγn− µ | − 1 = h2t γn| aγn− µ | + | aγn− µ | nγn − 1 ≥ h 2 t γn 1 (γn/2) + | aγn− µ | nγn − 1 by inequality (13) ≥ µ 2h2 t H2 min − 1 ¶ + | aγn− µ | nγn by inequality (20) > 0.
Therefore, it suffices to show that · nh2 t + (aγn− µ)2 nγn| aγn− µ | − 1 ¸2 − " p nh2 t + (aγn− µ)2 nγn #2 > 0.
Indeed, · nh2 t + (aγn− µ)2 nγn| aγn− µ | − 1 ¸2 − " p nh2 t + (aγn− µ)2 nγn #2 = 1 n2γ2 n " µ nh2 t | aγn− µ | + | aγn− µ | − nγn ¶2 − nh2 t − (aγn− µ)2 # = 1 n2γ2 n · n2h4 t − 2n2h2tγn| aγn− µ | (aγn− µ)2 + n2γ2 n− 2nγn| aγn− µ | + nh2t ¸ ≥ 1 n2γ2 n · n2h2 t (aγn− µ)2 (h2 t − γn2) + (n2− n) γn2+ nh2t ¸ by inequality (13) > 0, as desired.
Appendix C
We shall prove that the size of the MT tree is only quadratic in maturity N if criterion (23) holds. From Eq. (21),
η ≤ &p nh2 t + (γn2/4) nγn ' by inequality (13) = & 2pnh2 t√+ [ Hmin2 /(16n) ] n Hmin ' by equality (22) = &s 4h2 t H2 min + 1 4n2 ' ≤ » 2ht Hmin ¼ + 1 ≤ 2ht Hmin + 2. (26) Define C ≡ Hmin[ √ n + 1/(4√n) ] and D ≡ pβ1+ β2( √ n + c)2. According to
updating rule (15), h2t+1 ≤ max −n≤`≤n ( β0+ β1h2t + β2h2t · `ηγn+ aγn− (r − h2t/2) ht − c ¸2) = max −n≤`≤n n β0+ β1h2t + β2 £ `ηγn+ aγn− (r − h2t/2) − cht ¤2o ≤ max −n≤`≤n n β0+ β1h2t + β2 £ | `ηγn| + | aγn− (r − h2t/2) | + | − cht| ¤2o ≤ β0+ β1h2t + β2 ³ nηγn+ γn 2 + cht ´2 by inequality (13) ≤ β0+ β1h2t + β2 ·µ 2nht Hmin + 2n + 1 2 ¶ γn+ cht ¸2 by inequality (26) = β0+ β1h2t + β2 £ (√n + c) ht+ C ¤2 by equality (22) = £β1+ β2( √ n + c)2¤h2 t + 2β2C( √ n + c) ht+ β0+ β2C2 = (Dht+ P)2+ Q,
where P and Q are positive numbers independent of ht and t. Hence there exists
a number R > 0 independent of ht and t such that
h2
t+1 ≤ (Dht+ P)2+ Q ≤ (Dht+ R)2.
(For example, one can pick R = P +√Q.) This leads to
ht+1≤ Dht+ R. (27)
Let Ht stand for the largest of all the volatilities at date t. Note that H0 = h0.
By inequality (27), state (yt, h2t)’s 2n + 1 successor volatilities are all dominated by
Dht+R. Suppose Ht+1 at date t+1 follows volatility ˜ht at date t (which may be an
interpolated volatility). Then Ht+1≤ D˜ht+ R ≤ DHt+ R; hence Ht+1 ≤ DHt+ R.
By induction, Ht+1≤ R t X i=0 Di+ Dt+1h0 = R 1 − D + · h0+ R D − 1 ¸ Dt+1.
Therefore, Ht+1 at date t + 1 does not grow exponentially when D ≤ 1, i.e., when
β1+ β2( √
n + c)2 ≤ 1, which establishes criterion (23).
We proceed to derive the quadratic tree size. Assume D ≤ 1. It is important to note that Ht+1≤ R 1 − D + h0D t+1≤ R 1 − D + h0.
By bound (14) and inequality (26), the number of nodes at date t > 0 is at most t−1 X i=0 · 2n µ 2Hi Hmin + 2 ¶ + 1 ¸ = (4n + 1) t + 4n Hmin t−1 X i=0 Hi ≤ (4n + 1) t + 4n Hmin t µ R 1 − D + h0 ¶ = · 4n + 1 + 4n Hmin µ R 1 − D + h0 ¶ ¸ t.
Hence the number of nodes at date t is linear in t as depicted in Fig. 17. The total size of an N-day MT tree is therefore at most
1 + N X t=1 · 4n + 1 + 4n Hmin µ R 1 − D + h0 ¶ ¸ t = O(N2), a quadratic growth.
References
[1] Bali, Turan G. “An Empirical Comparison of Continuous Time Models of the Short Term Interest Rate.” The Journal of Futures Markets, 19, No. 7 (1999), 777–797.
[2] Bollerslev, T. (1986) Generalized Autoregressive Conditional Heteroskedas-ticity. Journal of Econometrics, 31, pp. 307–327.
[3] Bollerslev, T., Chou, R.Y., and Kroner, K. (1992) ARCH Modeling in Finance: a Review of the Theory and Empirical Evidence. Journal of
Econo-metrics, 52, pp. 5–59.
[4] Cakici, N., and Topyan, K. (2000) The GARCH Option Pricing Model: a Lattice Approach. Journal of Computational Finance, 3(4), pp. 71–85.
[5] Chien, H.-H. (2003) On the Complexity of the Ritchken-Sankarasubramanian
Interest Rate Model. MBA Thesis. Department of Finance, National Taiwan
University, Taiwan.
[6] Cox, J.C., Ingersoll, J.E., and Ross, S.A. (1985) A Theory of the Term Structure of Interest Rates. Econometrica, 53(2), pp. 385–407.
[7] Cox, J.C., Ross, S.A., and Rubinstein, M. (1979) Option Pricing: a Simplified Approach. Journal of Financial Economics, 7(3), pp. 229–263.
[8] Dai, T.-S. (2004) Pricing Asian Options with Lattices. Ph.D. Thesis. Depart-ment of Computer Science and Information Engineering, National Taiwan Uni-versity, Taiwan.
[9] Dai, T.-S., and Lyuu, Y.-D. (2004) An Exact Subexponential-Time Lattice Algorithm for Asian Options. In Proceedings of ACM-SIAM Symposium on
Dis-crete Algorithms, pp. 710–717, New Orleans, January 11–13, 2004. Philadelphia:
Society for Industrial and Applied Mathematics, 2004.
[10] Dempster, M.A.H., and Richards, D.G. (2000) Pricing American Options Fitting the Smile. Mathematical Finance, 10(4), pp. 157–177.
[11] Duan, J.-C. (1995) The GARCH Option Pricing Model. Mathematical
Fi-nance, 5(1), pp. 13–32.
[12] Duan, J.-C. (1997) Augmented GARCH(p, q) Process and Its Diffusion Limit.
Journal of Econometrics, 79, pp. 97–127.
[13] Duan, J.-C., Gauthier, G., Sasseville, C., and Simonato, J.-G. (2003) Approximating American Option Prices in the GARCH Framework. Journal of
Futures Markets, 23(10), 915–928.
[14] Duan, J.-C., and Simonato, J.-G. (2001) American Option Pricing under GARCH by a Markov Chain Approximation. Journal of Economic Dynamics
& Control, 22, pp. 1689–1718.
[15] Engle, R.F., and Ng, V. (1993) Measuring and Testing of the Impact of News on Volatility. Journal of Finance, 48, pp. 1749–1778.
[16] Engle, R.F., and Patton, A.J. (2001) What Good Is a Volatility Model?
Quantitative Finance, 1(2), pp. 237–245.
[17] Forsyth, P.A., Vetzal, K.R., and Zvan, R. (2002) Convergence of Nu-merical Methods for Valuing Path-Dependent Options Using Interpolation.
Re-view of Derivatives Research, 5, pp. 273–314.
[18] Heston, S.L., and Nandi, S. (2000) A Closed-Form GARCH Option Valu-ation Model. Review of Financial Studies, 13(3), pp. 585–625.
[19] Hull, J.C., and White, A. (1993a) Efficient Procedures for Valuing Euro-pean and American Path-Dependent Options. Journal of Derivatives, 1(1), pp. 21–31.
