<< стр. 12(всего 14)СОДЕРЖАНИЕ >>
5
0

1992 1993 1994 1995 1996 1997
X

3 month interest rates
3
2
Y*E-2
1
0

1992 1993 1994 1995 1996 1997
X

Figure 15.4. Interest rates time series: German (thick dotted line),
Japanese (thin dotted line), American (thick straight line), Thai (thin

self-п¬Ѓnancing strategy which produces a positive expected payoп¬Ђ:

вЂў at time t
вЂ“ borrow the portfolio (1+rj )в€’1 E(Оёj,t+h |Ft )Xj,t from the countries
вЂ“ lend (1 + r0 )в€’1 Yt to Thailand,
вЂ“ invest the diп¬Ђerence (1 + rj )в€’1 E(Оёj,t+h |Ft )Xj,t в€’ (1 + r0 )в€’1 Yt in
the numeraire currency at the risk-free rate r1 ,
вЂў at time t + h
вЂ“ withdraw the amount Yt+h from Thailand,
15.2 Estimating the coeп¬ѓcients of an exchange rate basket 335

E(Оёj,t+h |Ft )Xj,t+h ,
вЂ“ pay back the loan of
вЂ“ keep the diп¬Ђerence.

The expression for the proп¬Ѓt and for its expected value are:
p
= Yt+h в€’ E(Оёj,t+h |Ft )Xj,t+h
О t+h
j=1
пЈ« пЈ¶
p
(1 + rj )в€’1 E(Оёj,t+h |Ft )Xj,t в€’ (1 + r0 )в€’1 Yt пЈё
+(1 + r1 ) пЈ­
j=1
пЈ« пЈ¶
p
(1 + rj )в€’1 E(Оёj,t+h |Ft )Xj,t в€’ (1 + r0 )в€’1 Yt пЈё .
E(О t+h |Ft ) = (1 + r1 ) пЈ­
j=1

15.2.2 Estimation results

For the implementation of the investment strategy described above one needs
the estimate of the, possibly time-varying, basket weights. The precision of
the estimation has a direct impact on the economic result of the investment.
Therefore, we compare four diп¬Ђerent estimators of the basket weights: the
adaptive, the recursive, the window and the Kalman п¬Ѓlter estimator using
economic criteria for a one month and for a three month investment horizon.
In particular we compute the average expected proп¬Ѓt and the average realized
proп¬Ѓt.
The adaptive estimation procedure requires three parameters: m, О» and Вµ. The
choice of m0 does not inп¬‚uence the results very much and it can be reasonably
set to 30. This value represents the minimal amount of data which are used
for the estimation, and in the case of a structural break, the minimal delay
before having the chance of detecting the change point. The selection of О» and
Вµ is more critical. These two values determine the sensitivity of the algorithm.
Small values would imply a fast reaction to changes in the regressor coeп¬ѓcients,
but but they would also lead to the selection of intervals of homogeneity which
are possibly too small. Large values would imply a slower reaction and con-
sequently the selection of intervals which can be too large. To overcome this
problem we suggest the following approach.
The main idea is that small changes in the values of О» and Вµ should not aп¬Ђect
the estimation results. Therefore we restrict our attention on a set S of possible
336 15 Locally time homogeneous time series modeling

USD weight

8
0.028+Y*E-2
6
4
2

1994 1995 1996 1997
X

DEM weight
6
4
Y*E-2
2
0

1994 1995 1996 1997
X

JPY weight
0.7
0.6
0.5
Y
0.4
0.3
0.2

1994 1995 1996 1997
X
(straight line), recursive (thine dotted line), window (thick dotted line).

pairs (О», Вµ). In the present context we chose all the even number between 2
and 8:
S = {(О», Вµ)| О», Вµ в€€ {2, 4, 6, 8}}
15.2 Estimating the coeп¬ѓcients of an exchange rate basket 337

Then we compare the 16 pairs with the following criterion at each time t:
пЈ¶2
пЈ«
tв€’1 d
(О»в€— , Вµв€— ) = arg min пЈ­Ys в€’ Оёj,s|sв€’h Xj,s пЈё .
(О»,Вµ)в€€S s=tв€’200 j=1

Finally, we estimate the value of Оёt+h|t with the selected pair (О»в€— , Вµв€— ). The
appeal of the above selection criterion consists of the fact that it leads to the
choice of the pair (О», Вµ) which has provided the least quadratic hedging costs
over the past trading periods. Notice that in general we have diп¬Ђerent results
depending on the length of the forecasting horizon: here one and three month.
Figure 15.5 shows the results for the three month horizon. It is interesting to
see that the adaptive estimate tends to coincide with the recursive estimate
during the п¬Ѓrst half of the sample, more or less, while during the second half
of the sample it tends to follow the rolling estimate.
We remark that the problem of selecting free parameters is not speciп¬Ѓc to the
adaptive estimator. The window estimator requires the choice of the length
of the window: k, while the Kalman п¬Ѓlter needs the speciп¬Ѓcation of the data
generating process of Оёt and the determination of ОЈ and Пѓ. In this application
k is set equal to 250, ОЈ and Пѓ are estimated recursively from the data using the
OLS, while Оё0|0 and P0|0 are initialized using the п¬Ѓrst 350 observations which
are then discarded. We remark that this choice is consistent with the one of
Christoп¬Ђersen and Giorgianni (2000).
Table 15.2 shows the result of the simulated investment. The investments are
normalized such that at each trading day we take a short position of 100 USD
in the optimal portfolio of the hard currencies. The result refers to the period
April 9 1993 to February 12 1997 for the one month horizon investment and
June 7 1993 to February 12 1997 for the three month horizon investment. Notice
п¬Ѓrst that the average realized proп¬Ѓts are positive and, as far as the three month
investment horizon is concerned, they are signiп¬Ѓcantly larger than zero among
all methods. This provides a clear evidence for the fact that arbitrage proп¬Ѓts
were possible with in the framework of the Thai Bath basket for the period
under study. The comparison of the estimator also show the importance of
properly accounting for the time variability of the parameters. The recursive
estimator shows modest result as far as the realized proп¬Ѓts are concerned and
the largest bias between expected the realized proп¬Ѓt. On one side, the bias is
reduced by the window estimator and by the Kalman п¬Ѓlter, but on the other
side these two methods provide a worse performance as far as the realized proп¬Ѓt
are concerned. Finally, the adaptive estimator appears to be the best one, its
338 15 Locally time homogeneous time series modeling

ONE MONTH HORIZON Recursive Window KF Adaptive

Average Expected Proп¬Ѓts .772 .565 .505 .553

Average Realized Proп¬Ѓt .403 .401 .389 .420

Standard errors (.305) (.305) (.330) (.333)

THREE MONTH HORIZON Recursive Window KF Adaptive

Average Expected Proп¬Ѓts 1.627 1.467 1.375 1.455

Average Realized Proп¬Ѓt 1.166 1.141 1.147 1.182

Standard errors (.464) (.513) (.475) (.438)

Table 15.2. Summary statistics of the proп¬Ѓts.

bias is much smaller than the one of the recursive estimator and it delivers the
largest realized proп¬Ѓts for both investment horizons.

15.3 Estimating the volatility of п¬Ѓnancial time
series
The locally time homogeneous approach appears to be also appropriate for the
estimation of the volatility of п¬Ѓnancial time series. In order to provide some
motivation we п¬Ѓrst describe the stylized facts of п¬Ѓnancial time series. Let St
deп¬Ѓne the price process of a п¬Ѓnancial asset such as stocks or exchange rates,
then the returns are deп¬Ѓned as follows:

Rt = ln St в€’ ln Stв€’1 .

Stylized facts of п¬Ѓnancial asset returns are: a leptokurtic density, variance
clustering and highly persistent autocorrelation of square and absolute returns
(see Figure 15.6). Further details and examples on this topic can be found in
Taylor (1986) and in Franke, HВЁrdle and Hafner (2001).
a
15.3 Estimating the volatility of п¬Ѓnancial time series 339

returns
5
0
Y*E-2
-5

0 50 100 150 200 250 300
X

ACF of the absolute returns
15
10
Y*E-2
5
0

0 10 20 30 40 50
X
Figure 15.6. JPY/USD returns XFGretacf.xpl

15.3.1 The standard approach

The returns of п¬Ѓnancial time series are usually modeled by the following equa-
tion:
Rt = Пѓ t Оµt
Where Пѓt is a strictly positive process, which describes the dynamics of the vari-
ance of Rt , and Оѕt has a standard normal distribution: Оѕt в€ј N (0, 1). Standard
parametric models of the volatility are of (G)ARCH type:
2 2 2
Пѓt = П‰ + О±Rtв€’1 + ОІПѓtв€’1 ,
340 15 Locally time homogeneous time series modeling

like in Engle (1995) and Bollerslev (1995), and of stochastic volatility type:
2 2
ln Пѓt = Оё0 + Оё1 ln Пѓtв€’1 + ОЅt ,

as described by Harvey, Ruiz and Shephard (1995). These models have been
expanded in order to incorporate other characteristics of the п¬Ѓnancial return
time series: TARCH, EGARCH and QARCH explicitly assume an asymmet-
ric reaction of the volatility process to the sign of the observed returns, while
IGARCH and FIGARCH model the long memory structure of the autocorre-
lations of the square returns.

15.3.2 The locally time homogeneous approach

A common feature to all the models which have been cited in the previous
section is that they completely describe the volatility process by a п¬Ѓnite set of
parameters. The availability of very large samples of п¬Ѓnancial data has given
the possibility of constructing models which display quite complicated param-
eterizations in order to explain all the observed stylized facts. Obviously those
models rely on the assumption that the parametric structure of the process
remains constant through the whole sample. This is a nontrivial and possi-
bly dangerous assumption in particular as far as forecasting is concerned as
pointed out in Clements and Hendry (1998). Furthermore checking for param-
eter instability becomes quite diп¬ѓcult if the model is nonlinear, and/or the
number of parameters is large. Whereby those characteristics of the returns
which are often explained by the long memory and (fractal) integrated nature
of the volatility process, could also depend on the parameters being time vary-
ing. We want to suggest an alternative approach which relies on a locally time
homogeneous parameterization, i.e. we assume that the volatility Пѓ follows a
jump process and is constant over some unknown interval of time homogeneity.
The adaptive algorithm, which has been presented in the previous sections,
also applies in this case; its aim consists in the data-driven estimation of the
interval of time homogeneity, after which the estimate of the volatility can be
simply obtained by local averaging.

15.3.3 Modeling volatility via power transformation

Let St be an observed asset process in discrete time, t = 1, 2, . . . , П„ and Rt
are the corresponding returns: Rt = log(St /Stв€’1 ) . We model this process via
15.3 Estimating the volatility of п¬Ѓnancial time series 341

the conditional heteroscedasticity assumption

Rt = Пѓt Оµt , (15.10)

where Оµt , t в‰Ґ 1 , is a sequence of independent standard Gaussian random
variables and Пѓt is the volatility process which is in general a predictable
random process, that is, Пѓt is measurable with respect to Ftв€’1 with Ftв€’1 =
Пѓ(R1 , . . . , Rtв€’1 ) .
The model equation (15.10) links the volatility Пѓt with the observations Rt
via the multiplicative errors Оµt . In order to apply the theory presented in
Section 15.1 we need a regression like model with additive errors. For this
reason we consider the power transformation, which leads to a regression with
additive noise and so that the noise is close to a Gaussian one, see Carroll and
Ruppert (1988). Due to (15.10) the random variable Rt is conditionally on
Ftв€’1 Gaussian and it holds
2 2
E Rt |Ftв€’1 = Пѓt .

Similarly, for every Оі > 0 ,
Оі Оі Оі
E Rt |Ftв€’1 = Пѓt E (|Оѕ|Оі |Ftв€’1 ) = CОі Пѓt ,
2
Оі 2
Оі 2Оі 2Оі 2
= Пѓt E (|Оѕ|Оі в€’ CОі ) = Пѓt DОі
в€’ CОі Пѓt |Ftв€’1
E Rt

where Оѕ denotes a standard Gaussian r.v., CОі = E|Оѕ|Оі and DОі = Var|Оѕ|Оі .
2

Therefore, the process |Rt |Оі allows for the representation
Оі Оі
|Rt |Оі = CОі Пѓt + DОі Пѓt О¶t , (15.11)

where О¶t is equal (|Оѕ|Оі в€’ CОі ) /DОі . A suitable choice of the value of Оі provides
that the distribution of
(|Оѕ|Оі в€’ CОі ) /DОі
is close to the normal. In particular the value of Оі = 0.5 appears to be almost
optimal, see Figure 15.7.

15.3.4 Adaptive estimation under local time-homogeneity

The assumption of local time homogeneity means that the function Пѓt is
constant within an interval I = [П„ в€’ m, П„ ] , and the process Rt follows the
342 15 Locally time homogeneous time series modeling

Normal and powertransformed densities

1
Y
0.5
0

-1 0 1 2
X
Figure 15.7. Normal and power transformed densities for Оі = 0.5.
XFGpowtrans.xpl

Оі
regression-like equation (15.11) with the constant trend ОёI = CОі ПѓI which can
be estimated by averaging over this interval I :
1
|Rt |Оі .
ОёI = (15.12)
|I|
tв€€I

By (15.11)

CОі DОі 1 sОі
Оі Оі
ОёI = Пѓt + Пѓt О¶t = Оёt + Оёt О¶t (15.13)
|I| |I| |I| |I|
tв€€I tв€€I tв€€I tв€€I

with sОі = DОі /CОі so that

1
EОёI = E Оёt , (15.14)
|I|
tв€€I
2
s2 s2
Оі Оі 2
E Оёt О¶t = E Оёt . (15.15)
|I|2 |I|2
tв€€I tв€€I

Deп¬Ѓne also
s2
Оі
2 2
vI =2 Оёt .
|I|
tв€€I
15.3 Estimating the volatility of п¬Ѓnancial time series 343

In view of (15.15) this value is called the conditional variance of ОёI . Under
local homogeneity it holds Оёt is constantly equal to ОёI for t в€€ I , and hence,

EОёI = ОёI
s2 ОёI
2
Оі
2
vI = Var ОёI = .
|I|

A probability bound analogous to the one in Section 15.1 holds also in this
2
case. Let the volatility coeп¬ѓcient Пѓt satisfy the condition b в‰¤ Пѓt в‰¤ bB with
some constants b > 0, B > 1 . Then there exists aОі > 0 such that it holds for
every О» в‰Ґ 0

О»2
в€љ
P |ОёI в€’ ОёП„ | > О»vI в‰¤ 4 eО»(1 + log B) exp в€’ . (15.16)
2aОі

The proof of the statement above and some related theoretical results can be
found in Mercurio and Spokoiny (2000).
For practical application one has to substitute the unknown conditional stan-
dard deviation with its estimate: vI = sОі ОёI |I|в€’1/2 . Under the assumption of
time homogeneity within an interval I = [П„ в€’ m, П„ ] equation (15.16) allows
to bound |ОёI в€’ ОёJ | by О»vI + ВµvJ for any J вЉ‚ I, provided that О» and Вµ
are suп¬ѓciently large. Therefore we can apply the same algorithm described in
Section 15.1 in order to estimate the largest interval of time homogeneity and
the related value of ОёП„ . Here, as in the previous section, we are faced with
the choice of three tuning parameters: m0 , О», and Вµ. Simulation studies and
repeated trying on real data by Mercurio and Spokoiny (2000) have shown that
the choice of m0 is not particularly critical and it can be selected between 10
and 50 without aп¬Ђecting the overall results of the procedure.
As described in Section 15.2.2, the choice of О» and Вµ is more delicate. The
inп¬‚uence of О» and Вµ is similar to the one of the smoothing parameters in
the nonparametric regression. The likelihood of rejecting a time homogeneous
interval decreases with increasing О» and/or Вµ. This is clear from equation
(15.6). Therefore if О» and Вµ are too large this would make the algorithm too
conservative, increasing the bias of the estimator, while too small values of
О» and Вµ would lead to a frequent rejection and to a high variability of the
estimate. Once again, a way of choosing the optimal values of О» and Вµ can be
made through the minimization of the squared forecast error. One has to deп¬Ѓne
a п¬Ѓnite set S of the admissible pair of О» and Вµ. Then for each pair belonging
344 15 Locally time homogeneous time series modeling

(О»,Вµ)
to S one can compute the corresponding estimate: Оёt and then select the
optimal pair and the corresponding estimate by the following criterion:
T 2
(О»,Вµ)
|Rt |Оі в€’ Оёt
(О», Вµ) = min .
О»,Вµв€€S
t=0

Figure 15.8 shows the result of the on-line estimation of the locally time homo-
geneous volatility model for the JPY/USD exchange rate. The bottom plot, in
particular, shows the estimated length of the interval of time homogeneity: m,
at each time point.

15.4 Technical appendix
In this section we give the precise conditions under which the bound (15.4)
holds. Deп¬Ѓne:
WI = VIв€’1 ,
VI = Пѓ в€’2 Xt Xt
tв€€I

furthermore let wij,I denote theв€љ elements of WI . For some positive constants
b > 0, B > 1, ПЃ < 1, r в‰Ґ 1, О» > 2 and for i = 1 . . . p consider the random set
were the following conditions are fulп¬Ѓlled:
в€’1
пЈ± пЈј
b в‰¤ wjj,I в‰¤ bB
пЈґ пЈґ
пЈґ пЈґ
пЈґ пЈґ
пЈґ пЈґ
пЈґ пЈґ
пЈІ пЈЅ
||VI Вµ||2 wjj,I в‰¤ r
sup
Ai,I =
пЈґ {Вµв€€RK :||Вµ||=1} пЈґ
пЈґ пЈґ
пЈґ пЈґ
пЈґ пЈґ
пЈґ пЈґ
|wji,I /wjj,I | в‰¤ ПЃ в€Ђ i = 1, . . . , p
пЈі пЈѕ

Let (Y1 , X1 ) . . . (YП„ , XП„ ) obey (15.1), where the regressors are possibly stochas-
tic, then it holds holds for the estimate ОёI :
в€љ
P |Оёi,I в€’ Оёi,П„ | > О» wii,I ; Ai,I

в‰¤ 4e ln(4B)(1 + 2ПЃ r(d в€’ 1)О»)pв€’1 О» exp(в€’О»2 /2), i = 1, . . . , p.

A proof of this statement can be found in Liptser and Spokoiny (1999). For
a further generalization, where the hypothesis of local time homogeneity holds
only approximatively, see HВЁrdle et al. (2000).
a
15.4 Technical appendix 345

YEN/DM returns
10
5
Y*E-2
0
-5
-10

0 5 10 15 20
X*E3

Volatility
10
0.0002+Y*E-4
8
6
4
2

0 5 10 15 20
X*E3

Interval of homogeneity
10
Y*E2
5

0 5 10 15 20
X*E3
Figure 15.8. From the top: returns, estimated locally time homoge-
neous volatility and estimated length of the interval of time homogene-
ity. XFGlochom.xpl

Bibliography
Bollerslev, T. (1995). Generalised autoregressive conditional heteroskedasticity,
in Engle (1995).
346 15 Locally time homogeneous time series modeling

Carroll, R. and Ruppert, D. (1988). Transformation and Weighting in Regres-
sion, Chapman and Hall, New York.

Christoп¬Ђersen, P. and Giorgianni, L. (2000). Interest rate in currency basket:
Forecasting weights and measuring risk, Journal of Business and Eco-
nomic Statistics 18: 321вЂ“335.
Chui, C. and Chen, G. (1998). Kalman Filtering, Information Sciences, third
edn, Springer-Verlag, Berlin.
Clements, M. P. and Hendry, D. F. (1998). Forecastng Economic Time Series,
Cambridge University Press, Cambridge.
Cooley, T. F. and Prescott, E. C. (1973). An adaptive regression model, Inter-
national Economic Review 14: 364вЂ“371.
Eichengreen, B., Masson, P., Savastano, M. and Sharma, S. (1999). Transition
Strategies and Nominal Anchors on the Road to Greater Exchange Rate
Flexibility, number 213 in Essays in International Finance, Princeton Uni-
versity Press.

Elliot, R. J., Aggoun, L. and Moore, J. B. (1995). Hidden Markov Models,
Springer-Verlag, Berlin.
Engle, R. F. (ed.) (1995). ARCH, selected readings, Oxford University Press,
Oxford.
Franke, J., HВЁrdle, W. and Hafner, C. (2001). EinfВЁhrung in die Statistik der
a u
FinanzmВЁrkte, Springer, Berlin.
a

HВЁrdle, W., Herwartz, H. and Spokoiny, V. (2001). Time inhomogeous multiple
a
volatility modelling. Discussion Paper 7, Sonderforschungsbereich 373,
Humboldt-UniversitВЁt zu Berlin. To appear in Financial Econometrics.
a
HВЁrdle, W., Spokoiny, V. and Teyssi`re, G. (2000). Adaptive estimation for
a e
a time inhomogeneouse stochastic volatility model. Discussion Paper 6,
Sonderforschungsbereich 373, Humboldt-UniversitВЁt zu Berlin.
a
Harvey, A., Ruiz, E. and Shephard, N. (1995). Multivariate stochastic variance
models, in Engle (1995).
Lepski, O. (1990). One problem of adaptive estimation in gaussian white noise,
Theory Probab. Appl. 35: 459вЂ“470.
15.4 Technical appendix 347

Lepski, O. and Spokoiny, V. (1997). Optimal pointwise adaptive methods in
nonparametric estimation, Annals of Statistics 25: 2512вЂ“2546.

Liptser, R. and Spokoiny, V. (1999). Deviation probability bound for mar-
tingales with applications to statistical estimation, Stat. & Prob. Letter
46: 347вЂ“357.
Mercurio, D. and Spokoiny, V. (2000). Statistical inference for time-
inhomogeneous volatility models. Discussion Paper 583, Weierstrass In-
stitute for Applied Analysis and Stochastic, Berlin.
Mercurio, D. and Torricelli, C. (2001). Estimation and arbitrage opportunities
for exchange rate baskets. Discussion Paper 37, Sonderforschungsbereich
373, Humboldt-UniversitВЁt zu Berlin.
a
Musiela, M. and Rutkowski, M. (1997). Martingale Methods in Financial Mod-
elling, number 36 in Application of Mathemathics. Stochastic Modelling
and Applied Probability, Springer, New York.
Spokoiny, V. (1998). Estimation of a function with discontinuities via lo-
cal polynomial п¬Ѓt with an adaptive window choice, Annals of Statistics
26: 1356вЂ“1378.
Taylor, S. J. (1986). Modelling Financial Time Series, Wiley, Chichester.
16 Simulation based Option Pricing
Jens LВЁssem and JВЁrgen Schumacher
u u

16.1 Simulation techniques for option pricing
We introduce Monte Carlo techniques and Quasi Monte Carlo techniques for
option pricing. First, we give an idea how to use simulation techniques to
determine option prices, then - using the developed basic methods - we give
examples how to price more complex i.e. exotic options even on more than one
underlying. Finally we present a short guideline how to price exotic options
with the proposed techniques.
First, we take a look at a European put on one underlying stock, a pricing
problem which can be solved analytically e.g. by using the Black-Scholes for-
mula. We start with this problem not only because it has become a kind of
вЂќstandard problemвЂќ but also to have the possibility to compare the results of
our approximation with an analytical solution. At the same time we look at
the time-complexity of the used simulation technique. Next, we show how to
price path dependent options with Monte Carlo methods. Afterwards, we show
how to price a stock option on several underlyings. This implies that we have
to solve a multi-dimensional simulation problem.

16.1.1 Introduction to simulation techniques

The idea behind randomized algorithms is that a random sample from a pop-
ulation (of input variables) is representative for the whole population. As a
consequence, a randomized algorithm can be interpreted as a probability dis-
tribution on a set of deterministic algorithms.
We will see that there are three main advantages to randomized algorithms:
1. Performance: For many problems, it can be shown that randomized algo-
350 16 Simulation based Option Pricing

rithms run faster than the best known deterministic algorithm. 2. Simplicity:
Randomized algorithms are easier to describe and implement than comparable
deterministic algorithms. 3. Flexibility: Randomized algorithms can be easily
In general one distinguishes two types of randomized algorithms. Las Vegas
algorithms are randomized algorithms that always give correct results with
only the variation from one run to another being its running time. Monte
Carlo algorithms are randomized algorithms that may produce an incorrect
solution for which one can bound the probability of occurrence. The quality of
the solution can be seen as a random variable.
Within this chapter, we focus on Monte Carlo algorithms calculating the value
of the following integral
f (x)dx (16.1)
[0,1]d

by evaluation of f (x) for independent uniform distributed random vectors
X1 , X2 , . . . , Xn , Xi в€€ [0, 1]d .
The arithmetic mean of the values f (Xi ) can be seen as a guess for the expected
value of the random variable f (Xi ) and therefore can be interpreted as an
approximation for the value of the integral. According to the strong law of
large numbers the estimator for the expected value (the arithmetic mean of
the random function values) is converging to the expected value (the value of
the integral) with an increasing sample size. The probability that the absolute
error of the approximation result exceeds a п¬Ѓxed positive value is limited and
decreases to zero with an increasing sample size if the variance of f is п¬Ѓnite.

16.1.2 Pricing path independent European options on one
underlying

For the case of a European option on one underlying we have to approximate
the following integral via Monte Carlo simulation:

в€ћ
r(T в€’t)
CT (ST )g(ST |St , r, Пѓ, T в€’ t)dST (16.2)
e E [CT (ST )|St ] =
0

CT {f (x, St , r, Пѓ, T в€’ t)}dx
= (16.3)
[0,1)
16.1 Simulation techniques for option pricing 351

Where
2
)(T в€’t)))2
exp в€’ (log ST в€’(log St в€’(rв€’0.5Пѓ
2Пѓ 2 (T в€’t)
g(ST |St , r, Пѓ, T в€’ t) =
2ПЂПѓ 2 (T в€’ t)ST

is the risk neutral density function of the Black Scholes model with parameters:

ST : price of the underlying at maturity
St : price of the underlying at time t
r : risk free interest rate
Пѓ : volatility of log returns of the underlying
T в€’t : time to maturity

в€љ
1
ST = f (x, St , r, Пѓ, T в€’ t) = St exp{(r в€’ Пѓ 2 )(T в€’ t) + Пѓ T в€’ tF в€’1 (x)}
2
transforms the uniform distributed values x in g(ST |St , r, Пѓ, T в€’ t) distributed
underlying values ST . F в€’1 (x) is the inverse of the cumulative normal distri-
bution function and CT (y) is the payoп¬Ђ function of the option.
The Monte Carlo simulation calculates the value of the integral in the following
way:

1 n
1. n independent random underlying values ST . . . ST are generated by com-
puting f (x, St , r, Пѓ, T в€’t) for a set of uniformly distributed pseudo random
numbers X1 , . . . , Xn .
i i
2. the option payoп¬Ђ CT (ST ) is calculated for each ST .
3. the value of the integral in (16.3) is then approximated by the arithmetic
mean of the option payoп¬Ђs:
n
1
ВЇ i
C= CT (ST )
n i=1

We will now derive an estimate of the approximation error of the arithmetic
1 n
mean. We assume that ST . . . ST are independent random underlying samples
of the g(ST |St , r, Пѓ, T в€’ t) density. Using this assumption we can conclude that
ВЇ
C is a random variable with expected value
ВЇ = er(T в€’t) Ct (St )
E[C]
352 16 Simulation based Option Pricing

Additionally we have to assume that the variance of the option payoп¬Ђs CT (ST )
is given by:

2
CT (ST )2 g(ST |St , r, Пѓ, T в€’ t)dST в€’ E [CT (ST )]
Var [CT (ST )] = (16.4)
[0,в€ћ]

exists. Then we get:
n
1 1
ВЇ i
Var C = 2 Var CT (ST ) = Var [CT (ST )] (16.5)
n n
i=1

1 n
because of the independence of ST , . . . , ST .
ВЇ
The expected value of the random variable C equals the value of the inte-
gral er(T в€’t) Ct (St ) and its variance converges to zero with increasing n. The
probability that the approximation error is greater than a п¬Ѓxed positive value
decreases to 0 with an increasing number n. A п¬Ѓrst estimation of the error is
ВЇ
given by the Chebychev inequality for C,

1
n Var [CT (ST )]
ВЇ
P |C в€’ er(T в€’t) Ct (St )| в‰Ґ a в‰¤
a2
The bound given by this equation is rather imprecise since we do not make any
assumptions on the distribution of the random variable. Only the expected
value and the variance are used in the previous equation. According to the
ВЇ
central limit theorem the distribution of C converges to a normal distribution
for n в†’ в€ћ. It follows that the diп¬Ђerence between the approximation and the
ВЇ
integral, C в€’ er(T в€’t) Ct (St ) is approximately normally distributed with mean 0
and standard deviation

Var [CT (ST )]
ПѓC = (16.7)
ВЇ
n
for large n. According to Boyle (1977) a value of n > 1000 is suп¬ѓciently large
in order to use the normal distribution for error estimation purposes.

ВЇ
We get the following equation if we assume that C в€’ er(T в€’t) Ct (St ) is normal
distributed:

u2
1
ВЇ
C в€’ er(T в€’t) Ct (St ) в‰¤ =в€љ exp в€’
P du (16.8)
2ПѓC
2ПЂ ВЇ
в€’
16.1 Simulation techniques for option pricing 353

ВЇ
If we choose k as a multiple of the standard deviation ПѓC of C, then we get:
ВЇ

ВЇ
C в€’ er(T в€’t) Ct (St )
ВЇ
C в€’ er(T в€’t) Ct (St ) в‰¤ kПѓC в‰¤k
P =P
ВЇ
ПѓCВЇ
k
u2
1
=в€љ exp в€’ du
2
2ПЂ в€’k
=p (16.9)

в€љ
Given a п¬Ѓxed probability level p, the error converges to zero with O(1/ n).
The error interval holds for k = 1, 2, 3 with the respective probabilities p =
0.682, 0.955, 0.997
The conп¬Ѓdence intervals for a given probability level depend on the standard
deviation of the payoп¬Ђ function CT (ST ):

Пѓ CT = Var [CT (ST )] . (16.10)

In general, this standard deviation cannot be calculated with analytical meth-
ods. Therefore one calculates the empirical standard deviation Пѓ and uses it
ВЇ
as a proxy for the error bounds:

n
1 2
ВЇ
i
CT (ST ) в€’ C
Пѓ=
ВЇ . (16.11)
nв€’1
k=1

Figure 16.1 shows the evolution of the absolute error of the price for a European
call option calculated by Monte Carlo methods compared with в€љ analyticthe
solution. One can observe that the error tends to zero with O (1/ n).
We would like to give some of the main properties of algorithms using Monte
Carlo techniques. First from (16.9) it follows that the error bound tends to zero
в€љ
with O (1/ n) for a п¬Ѓxed probability level p. Second, the probability that a
в€љ
п¬Ѓxed error bound holds converges to 1 with O (1/ n), Mavin H. Kalos (1986).
Since these results hold independent of the dimension of the problem, which
aп¬Ђects only the variance of the payoп¬Ђ function with respect to the Black-Scholes
risk neutral density, the Monte Carlo method is especially well suited for the
evaluation of option prices in multidimensional settings. Competing pricing
methods e.g п¬Ѓnite diп¬Ђerences have exponential growing computational costs in
354 16 Simulation based Option Pricing

Errors in MC Simulation

0.03 0.06 0.09 0.12 0.15
absolute error

0 500000 1000000
number of iterations
Figure 16.1. Absolute error of a European Call option price calculated
by Monte Carlo simulations vs. nв€’1/2

the dimension of the problem. Another advantage of the Monte Carlo pricing
method is the error estimate given by the empirical standard deviation which
can be computed with a small additional eп¬Ђort.
The two most important drawbacks of the Monte Carlo simulation, mentioned
in literature are its small convergence speed compared to other techniques for
options on few underlyings and the diп¬ѓculties occurring for options with early
exercise possibilities. For example, American options giving the investor the
possibility to exercise the option at any time before and at maturity, are diп¬ѓcult
to price. To evaluate an American option means to п¬Ѓnd an optimal exercise
strategy which leads - using only basic Monte Carlo techniques - to a recursive
algorithm with exponential time-complexity. But more advanced techniques
using importance sampling methods show that Monte Carlo simulations can
be applied to evaluate American contracts, Broadie (2000).

16.1.3 Pricing path dependent European options on one
underlying

There are two categories of options. Path-independent options are options
whose payoп¬Ђ depends only on the underlying prices at maturity. Path-
16.1 Simulation techniques for option pricing 355

dependent options are options whose payoп¬Ђ depends on underlying price out-
comes St1 , . . ., Stm at several time points t1 в‰¤ . . . в‰¤ tm within the lifetime of
the respective option.
Within the group of path-dependent options one can distinguish options with
a payoп¬Ђ function depending on a continuously deп¬Ѓned path variable and op-
tions with a payoп¬Ђ function depending on a п¬Ѓxed number of underlying values.
The price of an option with many - usually equally spaced - exercise dates is
often approximated by the price of an option with a continuously deп¬Ѓned path
variable and vice versa.
Examples for path-dependent options are barrier options, lookback options,
and Asian options. The latter have a payoп¬Ђ function which is linked to the
average value of the underlying on a speciп¬Ѓc set of dates during the life of
the option. One distinguishes two basic forms of Asian options: options on
the geometric mean (for which the price can be calculated with standard tech-
niques) and options on the arithmetic mean (for which the price can not be
determined using standard approaches). Asian options are frequently used in
commodity markets. The volatility of the underlying prices of the commodities
is usually very high so that prices for vanilla options are more expensive than
for comparable Asian-style options.

16.1.4 Pricing options on multiple underlyings

In this section we show how to extend the Monte Carlo simulation technique
to higher dimensions. The problem is not only that one has to deal with higher
dimensional integrals, but also that one has to incorporate the underlying cor-
relation structure between the considered securities. In our framework we need
the covariance matrix of the log returns on an annual basis.
In general, a basket option is an option on several underlyings (a basket of
underlyings). Basket options can be European-, American or even Asian-style
options. Normally, the average of the underlying prices is taken to calculate
the price of the basket option, but sometimes other functions are used.
sional options is that the correlation between securities is taken into account.
This may lead to better portfolio hedges. We will look at a basket option on
п¬Ѓve underlyings where the underlying price of the best security in the basket
is taken to calculate the option price.
356 16 Simulation based Option Pricing

16.2 Quasi Monte Carlo (QMC) techniques for
option pricing

16.2.1 Introduction to Quasi Monte Carlo techniques

QMC methods can be considered as an alternative to Monte Carlo simulation.
Instead of (pseudo) random numbers, Quasi Monte Carlo algorithms use the
elements of low discrepancy sequences to simulate underlying values.
The discrepancy of a set of points P вЉ‚ [0, 1]s measures how evenly these points
are distributed in the unit cube. The general measure of discrepancy is given
by:

A(B; P )
в€’ О»s (B)
Dn (B; P ) := sup (16.12)
n
Bв€€B

where A(B; P ) is the number of points in P belonging to B, О»s (B) is the
Lebesgue measure of the set B, B is a family of Lebesgue measurable subsets
of [0, 1]s , and n is the number of elements in P .
The discrepancy of a set is the largest diп¬Ђerence between the number of points
in a subset and the measure of the subset. If we deп¬Ѓne B to be the family J of
s
subintervals i=1 [0, ui ), then we get a special measure, the star-discrepancy:
в€—
Dn (P ) := Dn (J ; P ) (16.13)

16.2.2 Error bounds

For the star-discrepancy measure and reasonable assumption on the nature of
the function that has to be integrated an upper bound on the error is given by
the following theorem:

THEOREM 16.1 (Koksma-Hlawka) If the function f is of п¬Ѓnite variation
V (f ) in the sense of Hardy and Krause, then the following equation holds for
all sets of points {x1 , . . . , xn } вЉ‚ I s = [0, 1]s
n
1 в€—
f (xi ) в€’ f (u)du в‰¤ V (f )Dn (x1 , . . . , xn ) (16.14)
n Is
i=1
16.2 Quasi Monte Carlo (QMC) techniques for option pricing 357

A proof is given in Niederreiter (1992).

This means that the error is bounded from above by the product of the
variation V (f ), which in our case is model and payoп¬Ђ dependent and the star-
discrepancy of the sequence. The bound cannot be used for an automatic error
estimation since the variation and the star-discrepancy cannot be computed
easily. It has been shown though that sequences exist with a star-discrepancy
of the order O(nв€’1 (ln n)s ). All sequences with this asymptotic upper bound
are called low-discrepancy sequences Niederreiter (1992). One particular
low-discrepancy sequence is the Halton sequence.

16.2.3 Construction of the Halton sequence

the interval [0, 1]. An element of this sequence is calculated by using the fol-
lowing equation:
в€ћ
nk,i pв€’kв€’1
xi = (16.15)
k=0
with i > 0, p = 2 and nk,i determined by the following equation:
в€ћ
nk,i pk ; 0 в‰¤ nk,i < p; nk,i в€€ N
i= (16.16)
k=0

Note that with the above equation nk,i is a function of i and takes values only
in {0; 1}. To illustrate the algorithm we calculate the п¬Ѓrst three points.

i = 1: n0,1 = 1, nk,1 = 0 for every k > 0
i = 2: n1,2 = 1, nk,2 = 0 for every k = 1
i = 3: n0,3 = n1,3 = 1, nk,3 = 0 for every k > 1

Therefore we get the sequence 1/2, 1/4, 3/4, 1/8, 5/8, .... The extension of this
construction scheme to higher dimensions is straightforward. For every dimen-
sion j = 1, . . . , d we deп¬Ѓne xj by
i
в€ћ
xj nk,i (j)pв€’kв€’1
= (16.17)
i j
k=0
358 16 Simulation based Option Pricing

with pj is the jth smallest prime number and nk,i (j) is calculated as follows:
в€ћ
nk,i (j)pk ; 0 в‰¤ nk,i (j) < pj ; nk,i (j) в€€ N в€Ђj
i= (16.18)
j
k=0

By using p1 = 2, p2 = 3 we get the following two-dimensional Halton sequence:
(1/2; 1/3), (1/4; 2/3), .... In contrast to grid discretization schemes like i/n i =
1, ..., n low-discrepancy sequences п¬Ѓll the integration space in an incremental
way avoiding the exponential growth of grid points of conventional schemes.
XploRe provides quantlets to generate pseudo random numbers and low dis-
crepancy sequences. For the generation of the pseudo random numbers we use

erg = randomnumbers (seqnum,d,n)
generates n pseudo random vectors of dimension d

where seqnum is the number of the random generator according to Table 16.1,
d is the dimension of the random vector and n the number of vectors generated.

0 Park and Miller with Bays-Durham shuп¬„e
1 LвЂ™Ecuyer with Bays-Durham shuп¬„e
2 Knuth
3 generator from G. Marsaglia et al. Marsaglia (1993)
4 random number generator of your system
5 generator from ACM TOMS 17:98-111
6 multiply with carry gen. (Marsaglia) Marsaglia (1993)
Table 16.1. Random generator that can be used in XploRe

The generation of low discrepancy sequences is provided by

erg = lowdiscrepancy (seqnum,d,n)
generates the п¬Ѓrst n low discrepancy sequence vectors of dimen-
sion d

where seqnum is the number of the low discrepancy sequence according to Table
16.2.
16.2 Quasi Monte Carlo (QMC) techniques for option pricing 359

0 Halton sequence
1 Sobol sequence
2 Faure sequence
3 Niederreiter sequence
Table 16.2. Low-discrepancy sequences available in XploRe,
(Niederreiter, 1992) .

16.2.4 Experimental results

Figure 16.2 shows that two dimensional Halton points are much more equally
spaced than pseudo random points. This leads to a smaller error at least for
вЂњsmoothвЂќ functions.

First 1000 random points First 1000 Halton points
1
1
dimension 2

dimension 2
0.5

0.5
0

0

0 0.5 1 0 0.5 1
dimension 1 dimension 1
Figure 16.2. 1000 two-dimensional pseudo random points
vs. 1000 Halton points XFGSOPRandomNumbers.xpl,
XFGSOPLowDiscrepancy.xpl

The positive eп¬Ђect of using more evenly spread points for the simulation task
is shown in Figure 16.3. The points of a low-discrepancy sequence are designed
in order to п¬Ѓll the space evenly without any restrictions on the independence
of sequence points where as the pseudo random points are designed to show no
statistically signiп¬Ѓcant deviation from the independence assumption. Because
of the construction of the low discrepancy sequences one cannot calculate an
360 16 Simulation based Option Pricing

Errors in QMC vs. MC Simulation

-5
log absolute error
-10

0 500000 1000000
number of iterations
Figure 16.3. Absolute error of a random sequence and the Halton
sequence for a put option

empirical standard deviation of the estimator like for Monte Carlo methods
and derive an error approximation for the estimation. One possible way out
of this dilemma is the randomization of the low-discrepancy sequences using
pseudo random numbers i.e. to shift the original quasi random numbers with
pseudo random numbers Tuп¬ѓn (1996). If x1 , . . . , xn are scalar elements of a
low-discrepancy sequence X then we can deп¬Ѓne a new low discrepancy sequence

xi + xi + <= 1
W ( ) = {y1 , . . . , yn } with yi = (16.19)
(xi + ) в€’ 1 otherwise

for a uniformly distributed value . Then we can calculate an empirical stan-
dard deviation of the price estimates for diп¬Ђerent sequences W ( ) for indepen-
dent values which can be used as a measure for the error. Experiments with
payoп¬Ђ functions for European options show that this randomization technique
reduces the convergence rate proportionally.
The advantage of the Quasi Monte Carlo simulation compared to the Monte
Carlo simulation vanishes if the dimension increases. Especially the compo-
nents with a high index number of the п¬Ѓrst elements in low-discrepancy se-
quences are not evenly distributed Niederreiter (1992). Figure 16.4 shows that
the 49th and 50th component of the п¬Ѓrst 1000 Halton points are not evenly
distributed. But the result for the п¬Ѓrst 10000 points of the sequence shows that
the points become more evenly spread if the number of points increases.
16.3 Pricing options with simulation techniques - a guideline 361

First 1000 Halton Points First 10000 Halton Points
1

1
dimension 50

dimension 50
0.5

0.5
0

0

0 0.5 1 0 0.5 1
dimension 49 dimension 49
Figure 16.4. The п¬Ѓrst 1000 and 10000 Halton points of dimension 49
and 50 XFGSOPLowDiscrepancy.xpl

However by using the Brownian Bridge path construction method we can limit
the eп¬Ђect of the high dimensional components on a simulated underlying path
and the corresponding path variable for the most common path dependent
options, Morokoп¬Ђ (1996). This method start with an empty path with known
start value and calculates at each step the underlying value for a time point
with maximum time distance to all other time points with known underlying
value until the whole path is computed. Experimental results show that we
can still get a faster convergence of the QMC simulation for options up to 50
time points if we apply this path construction method.

16.3 Pricing options with simulation techniques -
a guideline
In this section we would like to give a short guideline how to price exotic
options with Monte Carlo and Quasi Monte Carlo simulation techniques
within the framework described above. Furthermore we give some indications
about the limits of these techniques.
362 16 Simulation based Option Pricing

16.3.1 Construction of the payoп¬Ђ function

As a п¬Ѓrst step we have to deп¬Ѓne the payoп¬Ђ function corresponding to our
option product. Within the methods deп¬Ѓned in the quantlib finance we have
to consider three diп¬Ђerent cases.

One underlying + path independent
In this case the payoп¬Ђ function is called by the pricing routine with the sim-
ulated underlying value at maturity as the single argument. It calculates the
corresponding payoп¬Ђ and returns this value. We have deп¬Ѓned the payoп¬Ђ func-
tion for a put option with strike price 100 as an example for a one dimensional
payoп¬Ђ function.

Several underlying + path independent
For options whose payoп¬Ђ depends on the underlying values of several assets at
maturity, we have to deп¬Ѓne a payoп¬Ђ function on the vector of the underlying
values at maturity. An example for such an option is an exchange option that
permits to swap a deп¬Ѓned share with the best performing share in a basket. Its
payoп¬Ђ function is given by:

CT ((ST , . . . , ST )) = max{0, О±i (ST в€’ K i ) + 55 в€’ ST |i = 1, .., 5}
1 5 i 3

One underlying + path dependent
The third category of option types that are captured are path dependent op-
tions on one underlying. The payoп¬Ђ function of these options depends on the
underlying values at several п¬Ѓxed time points during the lifetime of the option.
Payoп¬Ђ functions for these contracts are called with a vector of underlying val-
ues whose ith element is the underlying value at the time ti which has to be
speciп¬Ѓed in the model.

16.3.2 Integration of the payoп¬Ђ function in the simulation
framework

After deп¬Ѓning the payoп¬Ђ function in XploRe we can start to calculate a price
estimate with the help of the appropriate simulation routine. In the one di-
mensional case we just have to call
16.3 Pricing options with simulation techniques - a guideline 363

erg = BlackScholesPathIndependent1D (s0,r,vola,dt,opt,
itr,gen)
MC estimation of the option price for a path independent option.
erg = BlackScholesPathIndependent1DQMC (s0,r,vola,dt,opt,
itr,gen)
QMC estimation of the option price for a path independent
option.

to get a price estimate and for the Monte Carlo case an empirical standard
deviation with respect to a start price of s0, a continuous risk free interest
rate of r, a volatility vola, a time to maturity of dt years, the payoп¬Ђ function
opt, sample size itr and the random/low-discrepancy generator with number
gen. Table 16.1 shows the random number generators and table 16.2 the low-
discrepancy generators that can be used. An application of these routines for
a Put option can be found in XFGSOP1DPut.xpl.
Pricing path-dependent options is only slightly more complicated. Here we
have to deп¬Ѓne the vector of time points for which underlying prices have to
be generated. This vector replaces the time to maturity used to price path
independent options. Then we can apply one of the following methods to
compute a price estimate for the path dependent option

erg = BlackScholesPathDependent1D (s0,r,vola,times,opt,
itr,gen)
MC estimation of the option price for path-dependent options.
erg = BlackScholesPathDependent1DQMC (s0,r,vola,times,opt,
itr,gen)
QMC estimation of the option price for path-dependent options,
with:

with respect to the start price s0, the continuous risk free interest rate r, the
volatility vola, the time scheme times, the payoп¬Ђ function opt, sample size
itr and the random/low-discrepancy generator with number gen, as given in
Tables 16.1 and 16.2. Using the above quantlets, we calculate the price of an
Asian call option in XFGSOP1DAsian.xpl.
364 16 Simulation based Option Pricing

In the case of multidimensional options we have to deп¬Ѓne a start price vector
and a covariance matrix instead of a single underlying price and volatility value.
Then we can call one of the multi-dimensional simulation routines:

erg = BlackScholesPathIndependentMD (s0,r,vola,dt,opt
,itr,gen)
MC estimation of the option price in the multidimensional Black
Scholes model

erg = BlackScholesPathIndependentMDQMC (s0,r,vola,dt,opt
,itr,gen)
QMC estimation of the option price in the multidimensional
Black Scholes model

with respect to the m dimensional start price vector s0, the continuous risk free
interest rate r, the mГ—m covariance matrix vola, the time to maturity dt, the
payoп¬Ђ function opt, the number of iterations itr and the generator number
gen according to the generators in Tables 16.1 and 16.2. Both quantlets are
illustrated in XFGSOPMD.xpl.
If in addition a dividend is paid during the time to maturity, we can use the
following two quantlets to calculate the option prices.

erg = BlackScholesPathIndependentMDDiv (s0,r,div,vola
,dt,opt,itr,gen)
MC estimation of the option price in the multidimensional Black
Scholes model
erg = BlackScholesPathIndependentMDDivQMC (s0,r,div,vola
,dt,opt,itr,gen)
QMC estimation of the option price in the multidimensional
Black Scholes model

The additional argument div is a m dimensional vector of the continuously paid
dividends. An application of these functions for our basket option is provided
in XFGSOPMDDiv.xpl.
 << стр. 12(всего 14)СОДЕРЖАНИЕ >>