University of Southern California Dissertations and Theses

Bayesian analysis of stochastic volatility models with Levy jumps

(USC Thesis Other)

Bayesian analysis of stochastic volatility models with Levy jumps

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content BAYESIAN ANALYSIS OF STOCHASTIC VOLATILITY MODELS WITH LEVY JUMPS by Pawel Szerszen A Dissertation Presented to the FACULTY OF THE GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (ECONOMICS) May 2008 Copyright 2008 Pawel Szerszen ii Acknowledgements I would never have been able to nish my dissertation without the guidance of my committee members and support from my wife and family. I would especially like to thank my Advisor, Professor Michael Magill and Pro- fessor Christopher Jones, for their generous time and commitment. Throughout my doctoral work they encouraged me to develop analytical thinking and research skills. Without their help this project would never have beeen accomplished. I am also very grateful for having an exceptional doctoral committee and wish to thank Roger Moon and Jaksa Cvitanic for their continual support and encour- agement. I am extremely grateful for their assistance, generosity, and advice. I also wish give special thanks to my wife Joanna for her support, understanding and patience during my PhD studies. Finally, I'd like to thank my parents. They were always supporting me and encouraging me with their best wishes. iii Table of Contents List of Tables iv List of Figures v Abstract vi Chapter 1. Introduction 1 1.1. L evy Processes 8 1.2. Variance Gamma Process 11 1.3. L evy stable Process 13 Chapter 2. Estimation Method 16 2.1. Model Specications 16 2.1.1. Dynamics of the Asset Returns Process 16 2.1.2. Discretization scheme 24 2.2. Markov Chain Monte Carlo methods 28 2.3. Bayesian Inference for Stable Distributions 32 2.4. Estimation Algorithm 36 2.5. Auxiliary Particle Filter 48 Chapter 3. Estimation Results 55 3.1. The Data 55 3.2. Identication problem 59 3.3. Estimation Results for S&P 500 Data 65 3.4. VaR analysis 80 Chapter 4. Conclusions 86 References 89 Appendix A. Derivation of complete conditional distributions 94 iv List of Tables Table 1. S&P500 daily log-returns descriptive statistics. 56 Table 2. Descriptive statistics of posterior distributions for models with in- nite activity jumps for simulated data. 64 Table 3. Descriptive statistics of posterior distributions for S&P500 return models. 65 Table 4. Posterior medians and 95% condence intervals. 76 Table 5. Descriptive statistics offb z (1) t g distribution. 81 Table 6. Descriptive statistics offb z (5) t g distribution. 81 Table 7. One-day ahead value at risk analysis (VaR) performance of models with stochastic volatility from diusion. 84 Table 8. Five-day ahead value at risk analysis (VaR) performance of models with stochastic volatility from diusion. 84 v List of Figures Figure 1. S&P 500 index log-level and log-returns. 57 Figure 2. Skewness and kurtosis of S&P500 returns for term structure of returns. 58 Figure 3. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (1) specication with stochastic volatility and no leverage eect, com- pound Poisson jumps. 67 Figure 4. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (2) specication with stochastic volatility and leverage eect, compound Poisson jumps. 68 Figure 5. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (3) specication with stochastic volatility from diusion, variance gamma jumps and no leverage eect. 69 Figure 6. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (4) specication with stochastic volatility from diusion, leverage eect and variance gamma jumps. 70 Figure 7. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (5) with stochastic volatility from variance gamma jumps. 71 Figure 8. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (6) specication with stochastic volatility from diusion, stable jumps and no leverage eect. 72 Figure 9. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (7) specication with stochastic volatility from diusion, leverage eect and stable jumps. 73 Figure 10. MCMC draws for S&P500 data (01/02/1981-02/02/2007) model (8) with stochastic volatility from stable jumps. 74 vi Abstract In this work we analyze asset returns models with diusion part and jumps in returns with stochastic volatility either from diusion or pure jump part. We consider dif- ferent specications for the pure jump part including compound Poisson, Variance Gamma and Levy -stable jumps. Monte Carlo Markov chain algorithm is con- structed to estimate models with latent Variance Gamma and Levy-stable jumps. Our construction corrects for separability problems in the state space of the MCMC for Levy-stable jumps. We apply our model specications for analysis of S&P500 daily returns. We nd, that models with innite activity jumps and stochastic volatility from diusion perform well in capturing S&P500 returns characteristics. Models with stochastic volatility from jumps cannot represent excess kurtosis and tails of returns distributions. One-day and one-week ahead prediction and VaR per- formance characterizing conditional returns distribution rejects Variance Gamma jumps in favor of Levy -stable jumps in returns. JEL classication: C1; C11; G1; G12 1 Chapter 1 Introduction Asset pricing theory and derivative pricing relies on the models of underlying se- curities. Foundations to the theory of stock prices or index levels modelling under both statistical and risk-neutral measures were set by Black-Scholes model in Black and Scholes (1973) ([6]). However, these models produce typically disappointing results as stated in Andersen, Benzoni and Lund (2001) ([2]). The Black-Scholes pricing formula based on the geometric Brownian motion model of stock prices produces biased estimates of derivative prices. Bates (1995) ([4]) argued, that the Black-Scholes formula tends to underprice the put options relative to call options. Other well documented biases are volatility smiles and skewness premia. There is an evidence reported by Rubinstein (1994) ([45]), that implied volatilities are not constant, as assumed in the Black-Scholes model, but vary across moneyness and maturity levels. The name "volatility smile" takes its roots in the specic form of dependence of implied volatilities on moneyness and maturity levels. The above ndings suggested further ways of asset pricing research including importance of jumps and stochastic volatility. Both of these elements were incor- 2 porated in the Andersen, Benzoni and Lund (2001) ([2]) paper with satisfactory results. Since derivative pricing is especially dependent upon volatility dynamics of the underlying security, the inclusion of stochastic volatility relaxes one of the main drawbacks of the Black-Scholes model. Moreover, as stated in Andersen, Benzoni and Lund (2001) ([2]), the extent of skewness, presence of outliers and excess kurto- sis in the underlying stock returns is also of great importance. Stochastic volatility and jumps, which help to induce skewness and outliers in the return distribution, are found to be good, but still partial solution of the problem. However, the gen- eralized version of the model required dierent techniques of estimation, where the special estimation problems come from the unobservable stochastic volatility pro- cess. This was nally resolved with development of ecient method of moments (EMM) procedure of Gallant and Tauchen (1996) ([21]) and Monte Carlo Markov Chain Methods. The next generalization called "leverage eect" allows for instantaneous corre- lation between increments of stock returns and increments of volatility. Empirical results of Jacquier, Polson and Rossi (2004) ([29]), Jones (2003) ([33]) and Ander- sen, Benzoni and Lund (2001) [2] found the respective correlation to be signicantly negative. The negative leverage eect has a deep intuitive explanation, since peri- ods of high volatility on the market coincides more often with market crashes. The leverage eect helps in capturing the skewness of the stock returns. Moreover, Jones 3 (2003) [33] included in one of his considered specications the leverage eect depen- dent upon level of volatility itself. His ndings suggest, that as volatility increases, the leverage eect is higher in magnitude. Hence, in the periods of high volatility the probability of market crashes is higher, than in periods with low volatility. One of the important questions concerns data used in the estimation of the model. There exists two general approaches to model the asset returns. First approach species the model under the statistical measure, which allows for direct analysis of the return series. The second one species the model under the risk- neutral measure and relies on the underlying options data. There is also a way to utilize information from both worlds as in Chernov and Ghysels (2000) [13], however, it results in even further technical diculties and has limited applicability in our time-changed framework. Moreover, estimation under both statistical and risk-neutral measures require denition of market premia for undiversiable risks as jumps and stochastic volatility. This causes joint estimation to involve additional possible sources of misspecication as noted by Andersen, Benzoni and Lund (2001) [2]. In this work we use data on daily S&P 500 index returns and perform estimation under statistical measure. We choose this data for empirical study, since this is a broad indicator of the equity market and has other comparable studies in the litera- ture. We consider model specications with dierent jumps in returns specications, 4 dierent source of stochastic volatility and either existence or lack of leverage eect. This enables us to study the marginal eects of jumps, source of stochastic volatility and leverage eect. Our choice for the jump structure includes standard compound Poisson jumps in returns and innite activity jumps with either nite or innite variation. In case of innite activity Levy jumps with nite variation we consider variance gamma process (VG) of Madan, Carr and Chang (1998) ([38]) and Levy stable jumps with index of stability 2 (1; 2) for innite activity jumps with innite variation. The source of stochastic volatility is taken to be either from diu- sion part or from the jump part in case of models with innite activity jumps. The choice of stochastic volatility from innite activity Levy jumps is, to our best knowl- edge, new to the literature with estimation under statistical measure and refers to time-changed methodology of Carr and Wu ([11]). Andersen, Benzoni and Lund (2001) ([2]) and Eraker, Johannes and Polson (2003) ([20]) shed a light on importance of Poisson jumps in the returns process. The discrete jumps in the returns cannot be modeled successfully by simple increase, or even jump, in the volatility from diusion. Since the arrival rate of Poisson jumps under statistical measure was found to be only around 1:5 jumps per year 1 , the more ne jumps cannot be modeled by infrequent and "large" in size compound Poisson jumps. This in turn is one of the main critiques of nite activity jumps in returns. 1 This estimate refers to model with square root volatility process. Please refer to Eraker, Johannes and Polson (2003) ([20]) for details. 5 The solution to this problem lies in introduction of innite activity Levy jumps in returns, that is the jumps with innite number of "small" jumps in a nite time interval. The latest generalizations include innite activity jumps as in the case of variance-gamma (VG) model of Carr et.al (2002) ([8]) and Li, Wells and Yu (2006) ([36]), which are also found to be important characteristics of the return process and perform better than Poisson type jumps. In this work we not only allow for innite activity jumps, but also for innite variation jumps by introduction of stable jumps in returns similarly to Huang and Wu (2004) ([27]) and Carr and Wu ([10]). In their models, however, the estimation is performed under risk-neutral measure and not under statistical measure. We ll this gap in the literature. Since stable jumps, except Gaussian case, do not have nite second moment, the CLT does not apply anymore. This gives the model ability to introduce non- normality on returns even for quarterly and yearly returns. This characteristic distinguishes Levy stable pure jump specication both from the nite activity Poisson jumps as in Andersen, Benzoni and Lund (2001) [2], and from the variance gamma pure jump process of Madan, Carr and Chang (1998) ([38]) and Li, Wells and Yu (2006) ([36]). With respect to the source of stochastic volatility our model specication draws from similar time-changed model of Huang and Wu (2004) ([27]), where the time- change have interpretation of business time versus calendar time. In their work, 6 however, authors estimate the parameters of the underlying S&P 500 returns models under risk-neutral measure with application of the fast Fourier transform (FFT) of Carr and Madan (1999) ([9]) used to nd the model implied option prices later to be compared with the observed options data. The time-change of the diusion, or pure jump part, is similar to introduction of stochastic volatility from diusion, or respectively pure jump part, as described in Carr and Wu ([11]). In this work we are able to address the problem of the source of stochastic volatility under statistical measure and answer the question if the continuous diusion or pure jump part is driven by the business or calendar time. We perform estimation under statistical measure using Monte Carlo Markov chain procedure. We face the problem of parameters estimation in presence of both latent volatility states and latent jump sizes. The recent attempt to estimate the model with latent stable jumps in returns in Li, Wells and Yu (2006) ([36]) intro- duces separability on the Markov chain state-space and hence leaves their results questionable. We ll this gap in the asset pricing literature by constructing the MCMC algorithm, which solves this problem. We also nd, that even using moder- ate sample size we can disentangle time-changed diusion from VG jumps and stable jumps under both stochastic volatility specications. The log-volatility specication in the stochastic volatility improves ability to estimate the parameters of the models with high degree of precision for the models with stable jumps. The models with 7 innite activity jumps in returns (both VG and stable jumps) are able to represent well excess kurtosis of the underlying data and tails of the returns distribution if diusion is the source of stochastic volatility. Models with stochastic volatility from pure jump part does not represent data well and this feature does not depend on the inclusion of the leverage eect from diusion, since models with stochastic volatility from diusion but no leverage eect perform signicantly better. Moreover, the leverage eect is found to be an important component of all models. In case of models with innite activity jumps it has an eect of higher degree of precision of parameter estimates by reducing the misspecication of independent increments of returns and log-volatility. In order to evaluate one-day and one-week ahead prediction performance of the models we construct an auxiliary particle lter algorithm for models with stochastic volatility from diusion, leverage eect and jumps. We nd, that models with VG jumps in returns perform the worst in terms of one-day and one-week ahead pre- diction, while models with Poisson and Levy stable jumps produce satisfactory results. This inability of VG jump models to represent one-day and one-week ahead conditional returns distribution can be explained by the lack of estimation precision of parameter governing skewness of returns. Since VG jumps are of innite activity, estimation precision in jump specic parameters is of special importance. Value at risk (VaR) analysis, which focuses on the left-tail characteristics of conditional 8 returns distribution, also conrms poor performance of models with VG jumps. Models with Poisson jumps and Levystable jumps perform satisfactory. In gen- eral, VG jump models tend to overestimate VaR values (in absolute value). Finally, the inclusion of leverage eect is found to increase VaR estimates (in absolute value) across all jump specications. 1.1 L evy Processes In this section we closely follow Applebaum (2004) [3] and Bertoin (1998) [5]. Let X = (X t ;t 0) be a scalar Levy process dened on a probability space ( ;z;P) with given ltration (z t ) t0 . From denition, a Levy process X has independent and stationary increments, or more precisely, X s X t is independent ofz t and has the same distribution asX st X 0 for all 0ts. It is also stochastically contin- uous. We restrict our analysis to the modication ofX which exhibits c adl ag paths and hence its sample paths are right-continuous with left limits 2 . By the Levy-Ito decomposition every Levy process can be decomposed as the sum of three indepen- dent processes: a linear drift, Brownian motion and a pure jump part. Accordingly, the log-characteristic function of a Levy process is a sum of log-characteristic func- tions of its Levy components and is given by the Levy-Khintchine formula. The 2 Please refer to Applebaum (2004) [3] for proof of existence of such c adl ag modication. 9 characteristic function of the Levy process is given by Xt (u) =E[exp(iuX t )] = exp(t x (u)) , t 0; (1.1) where u2 R and x (u) is called Levy or characteristic exponent of a given Levy process. Levy-Khintchine formula determines the functional form of this exponent: x (u) =ibu 1 2 u 2 2 + Z Rnf0g [exp(iuz) 1iuz1 jzj<1 ]v(dz); (1.2) where b2 R controls linear drift part, 0 controls the Brownian part and v is a Levy measure that characterizes the pure jump part of the Levy process. The triplet (b; 2 ;v) completely characterizes probabilistic behavior of the Levy process. The Levy measure v is a sigma-nite measure on Rnf0g; not necessarily a nite measure, satisfying Z Rnf0g min(1;z 2 )v(dz)<1: (1.3) The above condition implies nite quadratic variation of any Levy process. We can extend the Levy measurev to allR; without loss of generality, assumingv(f0g) = 0: The Levy measure has interpretation, that for any subset ER; v(E) is the rate at which Levy process takes jumps of size x2 E and measures numbers of jumps of size x2E in the unit time interval. The compound Poisson process is the only pure-jump ( = 0) Levy process which satisesv(R)<1, for any other pure jump processes v([";"]) =1 for any " > 0 and hence process exhibits innite number of small jumps in a nite time interval. For any Levy process, however, the number 10 of "large" jumps remains nite with v(E)<1 for any ER; E\f0g =;. From now on we call Levy process to have nite activity if and only if v(R) = <1. The value is then called Poisson intensity. This class of processes is very general and involves Brownian motion and com- pound Poisson process as two special cases, both already extensively studied in the asset returns literature. Brownian motion is the only Levy process with continu- ous sample paths and hence does not allow for discontinuous jumps. The compound Poisson jump process, however, represents special jump characteristics with its nite activity property. The sum of Brownian part and compound Poisson part although being a Levy Process, does not allow for more general jump structures and is one of the main reasons of critique of asset returns models based on them. In this work we allow for more general property of the jump structure by redening the jump part of the underlying asset returns process and allowing for innite activity. The pure jump Levy process with innite activity, however, can also be classied into two general sub-classes with respect to the total absolute variation of the process, that is the aggregate, absolute distance traveled by the process. The Levy pure jump process is of innite total variation if the following condition is satised by its Levy measure: 11 Z Rnf0g min(1;jzj)v(dz) =1; (1.4) otherwise it is called of nite total variation. The nite variation, innite activity jumps were already extensively studied in the literature with Madan, Carr and Chang (1998) ([38]) and Li, Wells and Yu (2006) ([36]) variance-gamma (VG) process. In this work, we concentrate on the comparison of dierent jump specications with respect to nite/innite activity and nite/innite total variation of the jump component. To represent the class of innite activity and nite variation jumps we choose the standard VG process of Madan, Carr and Chang (1998) ([38]). In case of innite variation, innite activity jumps we restrict our analysis to the Levy stable pure jump process with index of stability 2 (1; 2); which Levy measure satises condition (1.4) 3 . 1.2 Variance Gamma Process In this section we closely follow Madan, Carr and Chang (1998) ([38]). The VG pro- cess J VG is obtained by subordinating Brownian motion with drift j and variance 2 j at a time given by the independent gamma process G t with unit mean rate and variance rate j > 0. In the following let B t denote the standard Brownian motion 3 Please refer to Samorodnitsky and Taqqu (1994) [46] and Janicki and Weron (1994) [30] for further details. 12 on ( ;z;P), then J VG (t; j ; j ; j ) = j G t + j B Gt : (1.5) By specifying the triplet j 2 R; j > 0 and j > 0 we uniquely determine the probabilistic behavior of the process. As pointed out by Madan, Carr and Chang (1998) ([38]) the VG process gives the ability to control two characteristics of the distribution above that of volatility: the skewness is controlled by j and kurtosis by j . The characteristic function of the VG process is given by: E(exp(uiJ VG (t)) = 1 1i j j u + ( 2 j j =2)u 2 t= j : (1.6) The unconditional density of the VG process can be found by integrating out the time change componentG t . This complication, however, makes its direct application to the MCMC methods infeasible. The respective form of VG Levy measure relies on the fact, that VG process can be decomposed as the dierence of two independent gamma processes, that mean and variance rates are the functions of j , j and j . It can also be shown, that VG process is the process of innite activity and nite variation. For further details please refer to Madan, Carr and Chang (1998) ([38]) and Carr, Geman, Madan and Yor (2002) ([8]). 13 1.3 L evy stable Process In this section we closely follow Nolan (2005) ([41]) and Lombardi (2004) ([37]). The building block ofstable process is the stable distribution. LetS(;;; ) denotes stable distributed random variable with index of stability 2 (1; 2), skewness 2 [1; 1], scale parameter > 0, and location parameter2R. The main problem in application of stable distribution is the lack of density function except for few special cases, = 2 for normal distribution and = 1; = 0 for Cauchy distribution. It is not known, however, for index values 2 (1; 2). The family is identied by means of the characteristic function: E(exp(uiS) = 8 > > < > > : exp(iu juj [1isgn(u) tan 2 ]) for 6= 1 exp(iu juj[1 +i 2 sgn(u) log(u)]) for = 1: (1.7) The characteristic function (1.7) is the most widely used parametrization, for details see Samorodnitsky and Taqqu (1994) ([46]) and Zolotarev (1986) ([48]). If we set = 2 we have a normal distribution and with = 1 Cauchy distribution. In this parametrization controls the extent of skewness in the distribution with maximum positive skewness given by = 1 and maximum negative skewness given by =1. The extent of skewness disappears as % 2. There are also other alternative 14 specications, one of them given by the following characteristic function: E(exp(ui S) = 8 > > < > > : expfiu juj exp[i 2 sgn(u) min(; 2)]g for 6= 1 expfiu juj[1 +i 2 sgn(u) log( juj)]g for = 1; (1.8) where S(; ;; ) denote stable random variable with index of stability 2 (1; 2), skewness parameter 2 [1; 1], scale parameter > 0 and location parameter 2 R. It can be shown, that S(; ;; ) d = S(;;; ) if we use the following translation of the parameters: = cot( 2 ) tan( 2 min(; 2)) (1.9) = [cos( 2 min(; 2))] 1= ; while and remain unchanged. It is important to note, that for 2 (1; 2) the skewness parameter has a dierent interpretation, where > 0 denotes negative skewness and < 0 positive skewness. Moreover, the value of skewness is not a linear function of the skewness parameter and depends also on the index of stability. The parametrization (1.8), used in Buckle (1995) ([7]), is also applied in this work unless otherwise stated. It is important to note, that stable distribution is innitely divisible and hence there exists a Levy process with increments having stable distribution, the Levy stable motion 4 . 4 For proof and further details please see Applebaum (2004) ([3]). 15 The Levy stable motion, J SJ t , is a Levy process, where increments have cen- tered, stable distribution with unit scale: J SJ t J SJ s d =S(;; 0; (ts) 1= ) for 0st: (1.10) In this work we focus on the Levy stable process, which is obtained by pre- multiplying the Levy stable motion by a constant > 0; and using the scaling property 5 of the stable distribution, we have: (J SJ t J SJ s ) d =S(;; 0;(ts) 1= ) for 0st: (1.11) Since -scale parameter given by translation formula (1.9) is linear with respect to , there exists ; , such that following holds: (J SJ t J SJ s ) d =S(;; 0;(ts) 1= ) d = S(; ; 0; (ts) 1= ) for 0st; (1.12) so without loss of generality we can use both parametrizations if talking about increments of Levy stable process. Moreover, as pointed out in section (1.1), it can be shown, that Levy stable process J SJ t constructed above is a pure jump Levy process with innite activity and innite variation, for more detailed exposition please refer to Samorodnitsky and Taqqu (1994) ([46]) and Janicki and Weron (1994) ([30]). 5 Please see Nolan (2005) ([41]) for detailed discussion of stable distribution properties and their parametrizations. 16 Chapter 2 Estimation Method In this chapter we develop methodology to eciently estimate models with stochastic volatility, leverage eect and latent jumps in returns, the latent jumps with either Poisson, Variance Gamma or Levy stable specication. In the rst section we present all specications under consideration and their discretization schemes. In the next section we present foundations of Bayesian Econometrics and its Markov Chain Monte Carlo methods. Finally, we develop algorithm to estimate models with Levystable specication found to be the most challenging case. In the last section of this chapter we construct auxiliary particle lter algorithm for models with latent jumps in returns, stochastic volatility and leverage eect from diusion. 2.1 Model Specications 2.1.1 Dynamics of the Asset Returns Process Let Y t denote logarithm of asset price or logarithm of index level at time t and Y t+1 Y t be the corresponding log-return on the corresponding asset or index level. 17 We consider several specications, that dier by the source of stochastic volatility in the returns process and form of the jump component in returns. In the following (B (1) t ;B (2) t ) denes 2-dimensional standard Brownian motion on ( ;z;P) probabil- ity space dened above. Model (1) SVPJ: Stochastic Volatility from Diusion with compound Poisson jumps in returns. This model specication draws from Eraker, Johannes and Polson (2003) ([20]) with log-volatility specication as in Jacquier, Polson and Rossi (1994) ([28]) and is characterized by the following 2-dimensional stochastic dierential equation (SDE): dY t =dt + (exp(h t )) 0:5 dB (1) t +dJ PJ t (2.1) dh t = h ( h h t )dt + h dB (2) t ; where 2 R denes drift part of the return process, h 2 R denes the speed of mean reverting property of the log-volatility h t process towards its mean h 2 R and h > 0 denes the volatility of the log-volatility parameter. The compound Poisson jump process J PJ t is characterized by its Poisson arrival intensity j > 0 and normally distributed jumps with mean j 2R and variance 2 j : Note, that we do not include leverage eect and jumps in volatility in the SVPJ model. This specication allows for positive volatility of the underlying return process without introduction of some additional constraints on the parameters and is ex- 18 tensively used in the literature along with its extension for inclusion of jumps in volatility in Eraker, Johannes and Polson (2003) ([20]). In this and all other models with log-volatility process, we restrict the parameter values to ensure stationarity of the log-volatility process. As noted in Carr, Wu (2004) ([11]) stochastic volatility can be alternatively interpreted as the stochastic time change of the underlying processes, using the following specication: dY t =dt +dB (1) Tt +dJ PJ t (2.2) dh t = h ( h h t )dt + h dB (2) t ; T t = Z t 0 exp(h s )ds This characteristic gives the possibility to distinguish between the calendar time t and business timeT t driving the Brownian motionB (1) t , while the compound Poisson jumps are still driven by the calendar time. Model (2) SVPJLE: Stochastic Volatility from Diusion with Leverage Eect with compound Poisson jumps in returns. SVPJLE specication is a simple generalization of SVPJ model above with lever- age eect as in Jacquier, Polson and Rossi (2004) ([29]). Model SVPJLE is dened 19 by the 2-dimensional SDE: dY t =dt + (exp(h t )) 0:5 dB (1) t +dJ PJ t (2.3) dh t = h ( h h t )dt + h (dB (1) t + p 1 2 dB (2) t ); where2 [1; 1] denes instantaneous correlation between Brownian motions driv- ing returns process and log-volatility process. To illustrate this, let dene the B (3) t process: dB (3) t =dB (1) t + p 1 2 dB (2) t : (2.4) It can be shown, that B (3) t is a Brownian motion and E t (dB (1) t dB (3) t ) = dt and hence parameter 2 [1; 1] can be interpreted as the instantaneous correlation as noted above. The leverage eect has been extensively studied in the recent research and was found to be an important characteristic of the asset returns models both under sta- tistical measure (Jacquier, Polson and Rossi (2004) ([29]) and risk-neutral measure (Huang, Wu (2004) ([27])). Model (3) SVVG: Stochastic Volatility from Diusion with Variance Gamma jumps in returns. 20 This model specication draws from Li, Wells and Yu (2006) ([36]) and is dened by the following 2-dimensional SDE: dY t =dt + (exp(h t )) 0:5 dB (1) t +dJ VG t (2.5) dh t = h ( h h t )dt + h dB (2) t ; where J VG t is the VG process with parameters j 2 R, j > 0 and j > 0 as described in section (1.2). Note, that we do not allow for leverage eect in SVVG model. Model (4) SVVGLE: Stochastic Volatility from Diusion with Leverage Eect with Variance Gamma jumps in returns. Model SVVGLE is a simple generalization of model SVVG with inclusion of leverage eect as in Jacquier, Polson and Rossi (2004) ([29]) and is dened by the following SDE: dY t =dt + (exp(h t )) 0:5 dB (1) t +dJ VG t (2.6) dh t = h ( h h t )dt + h (dB (1) t + p 1 2 dB (2) t ); with all parameters dened as in models SVPJLE for the log-volatility equation and in model SVVG for the VG jump process. Model (5) VGSV: Stochastic volatility from VG Jumps only. 21 This model specication relies on Carr and Wu (2004) ([11]) time-changed Levy processes methodology. Model VGSV is dened by the following SDE, where the diusion part is driven by the calendar timet and the VG jump part by the business time T t : dY t =dt +dB (1) t +dJ VG Tt (2.7) d t = ( t )dt + dB (2) t T t = Z t 0 exp( s )ds where> 0 denes the volatility of diusion part,2R controls the deterministic drift,J VG t is the VG process and t is the "log-volatility" of the VG jump component. Similarly to the process for h t above is the mean of the log-volatility and controls the mean-reverting speed of the t process. Model (6) SVSJ: Stochastic Volatility from Diusion with Levy stable jumps in returns. In this model we replace the VG jumps in SVVG model and incorporate the innite activity, innite variation Levy stable jumps process. Let consider the following specication: dY t =dt + (exp(h t )) 0:5 dB (1) t + j dJ SJ t (2.8) dh t = h ( h h t )dt + h dB (2) t ; 22 which, using the Carr and Wu (2004) ([11]) and Huang and Wu (2004) ([27]), ac- counts for time changed diusion part and no time change in the pure jump of the process. The pure jump part of the process, given by j J SJ t ; ( j > 0), is a Levy stable jump process as dened in the previous section with given skewness pa- rameter 2 [1; 1] and index of stability 2 (1; 2). Please note, that we specify our model using second parametrization with increments distributed according to the characteristic function in (1.8). As pointed before, this is done without loss of generality, since we can always translate all the parameters back to the rst parametrization in (1.7). This model is similar to Li, Wells and Yu (2006) ([36]), however, we applied log- volatility specication and not level-volatility square-root model as in Heston (1993) ([26]). This specication helped to remedy several estimation issues. Moreover, the parameter governing the skewness is estimated and not xed as in Li, Wells and Yu (2006) ([36]). Model (7) SVSJLE: Stochastic Volatility from Diusion with Leverage Eect with Levy stable jumps in returns. 23 This specication extends the model (6) SVSJ by incorporating leverage eect and is given by the 2-dimesional system of SDE: dY t =dt + (exp(h t )) 0:5 dB (1) t + j dJ SJ t (2.9) dh t = h ( h h t )dt + h (dB (1) t + p 1 2 dB (2) t ): Model (8) SJSV: Stochastic volatility from Levy stable jumps only. This specication is the Black-Scholes model (1973) ([6]) with Levy stable pure jump part and stochastic volatility from jumps: dY t =dt +dB (1) t + exp( t )dJ SJ t (2.10) d t = ( t )dt + dB (2) t ; where model parameters are dened as in model (5) VGSV. In this model, similarly to model (5) VGSV, the only source of the stochastic volatility is the jump part of the returns process and not the diusion. The scaling property of the Levy stable process let us interpret the above model as the time-changed model, where the time change is applied only to the pure jump component J SJ t (Levy stable motion) and not to the diusion. Let dene the process jt : jt = exp( t ); (2.11) and hence: exp( t ) = 1= jt : 24 Then, using the Carr and Wu (2004) ([11]) time-change methodology, we have: dY t =dt +dB (1) t +dJ SJ Tt (2.12) d t = ( t )dt + dB (2) t T t = Z t 0 js ds = Z t 0 exp( s )ds where we use property (1.12) of Levystable process increments. This character- istic of the model SJSV gives us the possibility to distinguish between the calendar time t driving the diusion part and the business time T t driving the pure jump part J SJ t using simple specication in (2.10), where stochastic volatility from Levy stable jumps was introduced in the multiplicative way. Note, that this simple way of introducing stochastic volatility from jumps is not available in VGSV model. 2.1.2 Discretization scheme In order to estimate the parameters of the continuous-time specications we need to discreticize the models. In the following we use rst order Euler scheme 1 and " (1) t ;" (2) t ;" (J) t are independentiidN(0; 1) distributed, all other random variables are also independent: 1 Please refer to Kloeden and Platen (1992) ([35]) for details on higher order approximation techniques. 25 Model (1) SVPJ specication: Y t+ =Y t + + (exp(h t )) 0:5 " (1) t+ +q t+ t+ (2.13) h t+ =h t + h ( h h t ) + h p " (2) t+ ; q t+ iid Bernoulli( j ) and t+ iid N( j ; 2 j ) Model (2) SVPJLE specication: Y t+ =Y t + + (exp(h t )) 0:5 " (1) t+ +q t+ t+ (2.14) h t+ =h t + h ( h h t ) + h p (" (1) t+ + p 1 2 " (2) t+ ); q t+ iid Bernoulli( j ) and t+ iid N( j ; 2 j ) Model (3) SVVG specication: Y t+ =Y t + + (exp(h t )) 0:5 " (1) t+ +J VG t+ (2.15) h t+ =h t + h ( h h t ) + h p " (2) t+ ; J VG t+ = j G t+ + j p G t+ " (J) t+ ; G t+ iid gamma( j ; j ) 26 Model (4) SVVGLE specication: Y t+ =Y t + + (exp(h t )) 0:5 " (1) t+ +J VG t+ (2.16) h t+ =h t + h ( h h t ) + h p (" (1) t+ + p 1 2 " (2) t+ ); J VG t+ = j G t+ + j p G t+ " (J) t+ ; G t+ iid gamma( j ; j ) Model (5) VGSV specication: Y t+ =Y t + + p " (1) t+ +J VG t+ (2.17) t+ = t + ( t ) + p " (2) t+ J VG t+ = j G t+ + j p G t+ " (J) t+ G t+ jf t gindependent gamma( exp( t ) j ; j ): Model (6) SVSJ specication: Y t+ =Y t + + (exp(h t )) 0:5 " (1) t+ +S t+ (;; 0; j 1= ) (2.18) h t+ =h t + h ( h h t ) + h p " (2) t+ ; Model (7) SVSJLE specication: Y t+ =Y t + + (exp(h t )) 0:5 " (1) t+ +S t+ (;; 0; j 1= ) (2.19) h t+ =h t + h ( h h t ) + h p (" (1) t+ + p 1 2 " (2) t+ ); 27 Model (8) SJSV specication: Y t+ =Y t + + p " (1) t+ + exp( t )S t+ (;; 0; 1= ) (2.20) t+ = t + ( t ) + p " (2) t+ : In the above t = 1; 2:::;T; S t (;; 0; ) denotes i.i.d. centered stable random variables with index of stability2 (1; 2), skewness coecient2 [1; 1] and scale parameter > 0 in the parametrization given by the characteristic function in (1.8) as in Buckle (1995) ([7]). All other variables and parameters are dened in section (2.1.1). The problem we face concerns choice of> 0 parameter governing the extent of the discretization bias. In this work we x = 1 and use the data at daily frequency. As noted by Li, Wells and Yu (2006) ([36]) and Eraker, Johannes and Polson (2003) ([20]) the discretization bias of daily data should not be signicant, however, we leave other choices of < 1 for the further research, where the data points at not observed frequencies are treated as latent variables, for details please refer to Jones (2003) ([33]) and Jones (1998) ([32]). In models (6)-(8) we follow Buckle (1995) ([7]) specication of stable distribution characteristic function as in eq. (1.8), which for2 (1; 2) has negative skewness for > 0 and maximum negative skewness at = 1. This allows us to model negative skewness in asset returns. 28 2.2 Markov Chain Monte Carlo methods In this section we brie y describe Markov Chain Monte Carlo (MCMC) methods, more detailed exposition can be found in Chib and Greenberg (1996) ([14]), Johannes and Polson (2003) ([31]) and Jones (1998) ([32]). LetY =fY t g T t=1 denote the observed prices,X are the unobserved (latent) state variables and are the parameters of the model, whereX t = (X 1 ;X 2 ;:::;X L )2R L ; = ( 1 ; 2 ;:::; K ) 2 R K and L;K 2 N denote respectively dimensions of the state variable and parameter vectors. In the Bayesian inference we utilize the prior information on the parameters to derive the joint posterior distribution for both parameters and state variables. By the Bayes rule, we have: p(;XjY )/p(YjX;) p(Xj) p(); (2.21) where p(YjX;) is the likelihood function of the model, p(Xj) is the probability distribution of state variables conditional on the parameters and p() is the prior probability distribution on the parameters of the model. Ideally we would like to know the analytical properties of the joint posterior distribution ofX and, however, it is hardly feasible. The highly multidimensional joint posterior distribution is very often too complicated to work with and analytically intractable and hence even the direct simulation from the joint posterior distribution is also hard to perform. The remedy to this problem is the process of breaking the multidimensional distribution 29 p(;XjY ) into its complete conditional distributions proposed by Cliord and Ham- mersley, who proved, that the set of complete conditional distributions completely characterize the joint distribution. In other words, knowledge of conditional distribu- tions p(jX;Y ) and p(Xj;Y ) determine the joint posterior distribution p(;XjY ). We can continue in this manner and characterize the joint posterior distribution p(;XjY ) by the following set of complete 1-dimensional conditional distributions: p( 1 j 2 ; 3 ;:::; K ;X;Y ); p( 2 j 1 ; 3 ;:::; K ;X;Y ),..., p( K j 1 ; 2 ;:::; K1 ;X;Y ) and for the state variables p(X 1 jX 2 ;X 3 ;:::;X L ;;Y ), p(X 2 jX 1 ;X 3 ;:::;X L ;;Y ), . . . , p(X L jX 1 ;X 2 ;:::;X L1 ;;Y ). However, we can always group the variables in the several blocks if we have knowledge on the respective higher dimensional conditional distributions. The mild positivity conditions under which the Cliord and Hammer- sley theorem holds require, that for each point in the sample space, p(;XjY ) and the marginal distributions have strictly positive mass. The MCMC algorithm can be dened as a way to construct a Markov Chain, which invariant distribution is the desired target distribution. The simplest MCMC algorithm is a Gibbs sampler developed in Geman and Geman (1984) ([22]). The Gibbs sampler, under mild technical conditions, gives the possibility to construct draws from non-standard distributions using its complete conditional distributions. The Gibbs sampler is the following iterative procedure: given ( (n1) ;X (n1) ) at the 30 beginning of step n = 1;:::N Step n:1:1 Draw (n) 1 from p( 1 j (n1) 2 ; (n1) 3 ;:::; (n1) K ;X (n1) ;Y ) Step n:1:2 Draw (n) 2 from p( 2 j (n) 1 ; (n1) 3 ;:::; (n1) K ;X (n1) ;Y ) Step n:1:K Draw (n) K from p( 2 j (n) 1 ; (n) 2 ;:::; (n) K1 ;X (n1) ;Y ) Step n:2:1 Draw X (n) 1 from p(X 1 jX (n1) 2 ;X (n1) 3 ;:::;X (n1) L ; (n) ;Y ) Step n:2:2 Draw X (n) 2 from p(X 2 jX (n) 1 ;X (n1) 3 ;:::;X (n1) L ; (n) ;Y ) Step n:2:L Draw X (n) L from p(X L jX (n) 1 ;X (n) 2 ;:::;X (n) L1 ; (n) ;Y ) Continuing in this fashion, we construct sequence of draws ( (n) ;X (n) ) n=1;:::N , which converges to the stationary distribution of the Markov Chain given by the target densityp(;XjY ). As noted above, we can group several variables and draw from the respective conditional distributions. The proof of the Gibbs sampler and sucient conditions can be found in Chib and Greenberg (1996) ([14]). The mild sucient conditions refer to the irreducibility and aperiodicity of such constructed chain. It has to be guaranteed, that the algorithm can "visit" all values, that are in the support of the desired target density and that the support cannot be separated into several regions, that it could never leave. The separate problem is the speed of convergence, which in general cases is hard to be determined analytically, for 31 discussion on the speed of convergence please refer to Johannes and Polson (2003) ([31]) and Roberts and Rosenthal (1998) ([44]). In this work we rely on empirical verication if the draws produced by the algorithm, from some point on, can be treated as draws from the target, invariant distribution. The Gibbs sampler provides useful method to draw samples from complicated and non-standard distributions, however, it assumes, that we can sample directly from the set of all complete conditional distributions. If we face the problem of sampling from intractable distribution we can replace the particular Gibbs sampler step by the Metropolis-Hastings (MH) in Metropolis et al. (1953) ([40]). The gener- alized version of the MH algorithm can be summarized by the following procedure: given x (n1) at the beginning of step n = 1;:::N Step 1 Draw x (n) from the proposal density q(x (n) jx (n1) ) (2.22) Step 2 Accept x (n) with probability (x (n1) ;x (n) ); where (x (n1) ;x (n) ) = min ~ p(x (n) )q(x (n1) jx (n) ) ~ p(x (n1) )q(x (n) jx (n1) ) ; 1 : (2.23) In the above x can denote any state variable or parameter, whose complete condi- tional distribution ~ p in the Gibbs sampler step was found to be intractable. The speed of convergence of the algorithm depends heavily on the choice of the proposal density, from which we can sample without diculty. Under mild conditions on both the proposal density and target density ~ p the chain converges to its stationary 32 distribution given by the target density ~ p, please refer to Chib and Greenberg (1996) ([14]) for further details. If we apply the MH step within the Gibbs sampler, we do not have to iterate on MH-step in order to obtain draw from the desired ~ p distribu- tion, but rather perform single draw in each step of the Gibbs sampler. The Markov chain constructed in this way preserves the desired target density as its invariant distribution as proved in Chib and Greenberg (1996) ([14]). In our work we are interested in obtaining random samples from the posterior distribution p(;XjY ). This allows for computation of several statistics including the sample means and higher moments from the desired marginal posterior distri- butions. The sample mean from the posterior distribution of parameters is taken to be the population parameter estimate. Moreover, the ergodic averaging theorem guarantees almost sure convergence to the true population moments, for proof please see Johannes and Polson (2003) ([31]). 2.3 Bayesian Inference for Stable Distributions The main problem in application of the Bayesian MCMC methods for stable distri- butions is the nonexistence of its density function for index of stability 2 (1; 2). Buckle (1995) ([7]) found a solution to this problem by introducing auxiliary vari- able, such that the joint density of the auxiliary variable and the stable distributed random variable exists. Let ~ S and ~ Y be the random variables with their joint density 33 f, conditional on , , , ; f : (1; 0) ( 1 2 ;l ; )[ (0;1) (l ; ; 1 2 )! (0;1) given by: f(~ s; ~ y) = j 1j exp ( z t ; (~ y) =(1) ) z t ; (~ y) =(1) 1 jzj ; (2.24) where z = (~ s) ; t ; (~ y) = sin[~ y + ; ] cos~ y cos~ y cos[( 1)~ y + ; ] (1)= ; (2.25) and 2 (1; 2), 2 [1; 1], 2 (1; +1); > 0; ~ s2 (1; +1), ~ y2 ( 1 2 ; 1 2 ), with ; =(2)=2 and l ; = ; =. According to theorem 1. in Buckle (1995) ([7]) f is a proper joint density of ( ~ S; ~ Y ) and the marginal distribution of ~ S is S(;;; ). It is important to note, that the domain of the density function is (~ s; ~ y)2 (1; 0)( 1 2 ;l ; )[(0;1)(l ; ). Hence we have the following dependence between ~ S and ~ Y random variables: ~ S > 0 () ~ Y2 (l ; ; 1 2 ) (2.26) and ~ S < 0 () ~ Y2 ( 1 2 ;l ; ): (2.27) In our application we have to draw S t for all t conditional on all other state variables and parameters as in the Gibbs sampler. Since one of the conditioning state variables occurs to be the respective auxiliary variable ~ Y t ; it uniquely determines 34 the sign of the draw S (i) t at the ith step of the Gibbs sampler. This violates one of the main assumptions of the MCMC method, since the state space cannot be separated into two subspaces according to the starting value of the sign of S (1) t ; that it would never leave. To illustrate the problem, let the starting values in the Gibbs sampler specifyS (1) t > 0, ~ Y (1) t 2 (l ; ; 1 2 ) for somet, and all other parameters, including (1) , (1) , (1) , (1) (consistent with chosen values S (1) t and ~ Y (1) t , in the support of the joint distribution) and of other state variables. Suppose, without loss of generality, that we have to update the jump size S t as rst in the algorithm, then since ~ Y (1) t 2 (l ; ; 1 2 ) we have to draw S (2) t > 0, and later, all other jump specic parameters (2) , (2) , (2) , (2) consistent with S (2) t > 0; ~ Y (1) t 2 (l ; ; 1 2 ). At the end we have to update also the auxiliary variable ~ Y (2) t in the support of the joint distribution, hence ~ Y (2) t 2 (l ; ; 1 2 ). Continuing in this manner we constructed an MCMC chain, that never "visits" negative values of jump sizes at time t. The algorithm has to drawS (i) t for all iterationsi with the same sign as the starting value S (1) t . The above state-space separability problem is not an issue, however, if we do not treat the jump variables as the latent variables to be estimated and we observe the jump sizes S t for all t as in the Buckle (1995) ([7]) application simply, because we do not have to update observed jump sizes in the Gibbs sampler algorithm. Li, Wells and Yu (2006) ([36]) do not correct for the separability problem in their MCMC algorithm derived for the latent stable jumps. This leaves their results 35 questionable and demands alternative approach to the estimation of the latent stable distributed jumps. In this work we oer solution to this problem by construction of the mixture distribution of two, truncated at zero, stable distributions. Let dene the following probability: p ;; =P ( ~ S > 0; ~ Y2 (l ; ; 1 2 )) =P ( ~ S > 0); (2.28) which for = 0 can be found analytically to be: p ; =p ;;=0 = (2.29) 0:5 + arctan[ ^ tan(=2)] ; where ^ is given by ^ = cot( 2 ) tan( 2 min(; 2)) the translation of parameters formula in eq. (1.9) to Buckle (1995) ([7]) characteristic function specication. This formula is based on the value of the stable distribution function at the origin for characteristic function specication in eq. (1.7) as in Nolan (2005) ([41]).In the next step we have to dene the distribution of truncated stable variables ~ S + , ~ S and their respective auxiliary variables ~ Y + , ~ Y . In the following let ( ~ S + ; ~ Y + ) and ( ~ S ; ~ Y ) have joint distributions dened respectively by the following density functions: f + (~ s + ; ~ y + ) = 8 > > < > > : 1 p ; f(~ s + ; ~ y + ; = 0) if ~ s + > 0 and ~ y + 2 (l ; ; 1 2 ) 0 otherwise ; (2.30) f (~ s ; ~ y ) = 8 > > < > > : 1 1p ; f(~ s ; ~ y ; = 0)if ~ s < 0 and ~ y 2 ( 1 2 ;l ; ) 0 otherwise ; (2.31) 36 where the densityf is dened in eq. (2.24). Moreover, let ~ U be Bernoulli distributed (conditional on and ) random variable with probability of success p ; . Let ~ U, ( ~ S + ; ~ Y + ), ( ~ S ; ~ Y ) be independent, then it is straightforward to show, that: ~ S =S(;; = 0; ) d = ~ U ~ S + + (1 ~ U) ~ S : (2.32) This specication complicates the MCMC algorithm by introducing mixing variable ~ U and separate auxiliary variables ~ Y + and ~ Y for respectively positive jumps ~ S + and negative jumps ~ S , all of them to be updated in the MCMC algorithm. However, the above specication solves the problem of separability of the resulting Markov chain state-space in the latent stable jumps models of asset returns. 2.4 Estimation Algorithm In this section we brie y describe the set of complete conditional distributions to be used in the Gibbs sampler algorithm. Since models with Poisson jumps and VG jumps in returns with stochastic volatility from diusion were already studied in the literature, we postpone their model specic derivations to the appendix 2 . Note, that model VGSV in eq. (2.17) can also be treated in the similar way to models SVVG and SVVGLE by conditioning on the latent log-volatility from jumpsf t g. 2 Please refer to Li, Wells and Yu (2006) ([36]) for further details. 37 The update of latent state variablesf t g in the VGSV model is also postponed to the appendix. In the following we concentrate on the jump specic parameters and state vari- ables in the models with Levy stable jumps, that is the models SVSJ, SVSJLE and SJSV in equations (2.18), (2.19) and (2.20) respectively. In the sequel we assume the number of daily observations T and discretization parameter = 1. 1. Models SVSJ (2.18) and SVSJLE (2.19). Note, that the only dierence between SVSJ and SVSJLE specications is exis- tence of the leverage eect, hence it is enough to x = 0 for the SVSJ model. In the following, we present the detailed discussion of updating pure jump sizesf ~ S + t g T t=2 , f ~ S t g T t=2 , their respective auxiliary variablesf ~ Y + t g T t=2 ;f ~ Y t g T t=2 , the mixing variables f ~ U t g T t=2 and the jump specic parameters , ,, j . Let = nfg for sets and , S t = ~ U t ~ S + t + (1 ~ U t ) ~ S t and " (3) t = " (1) t + p 1 2 " (2) t . Let = (; h ; h ; h ;;;; j ) be the vector of parameters, X = (fh t g T t=1 ;f ~ U t g T t=2 ;f ~ S + t g T t=2 ;f ~ S t g T t=2 ;f ~ Y + t g T t=2 ;f ~ Y t g T t=2 ) (2.33) be the vector of state variables and 38 Y = (Y 1 ;Y 2 ;:::;Y T ) =fY t g T t=1 be the vector of observed log-asset prices or log-index levels. updating auxiliary variables for positive and negative part of stable jumps ~ Y + t for t = 2;:::;T Dene the following change of variables v + t = t ; ( ~ Y + t ). As proved in Buckle (1995) ([7]) function t ; : ( 1 2 ; 1 2 ) 7! R in eq. (2.25) is increasing for given parameters 2 (1; 2) and 2 [1; 1]. Moreover, for 2 (1; 1), t ; %1 as ~ Y + t % 1 2 andt ; &1 as ~ Y t & 1 2 . For the proof of these properties please refer to Buckle (1995) ([7]). Since we are interested in the negative skewness > 0, we have the following proposition for the maximum negative skewness = 1: Proposition 1. For 2 (1; 2) and = 1 maximum negative skewness, we have t ; ( ~ Y + t )% 1 1 (1)= as ~ Y + t % 1 2 . Proof. Apply L'Hospital's rule for cos~ y cos[(1)~ y+ ; ] for = 1: Then substitute the limit to formula for t ; in eq. (2.25). The result follows immediately. From (2.24) the conditional posterior for ~ Y + t is given by: p( ~ Y + t j;X ( ~ Y + t ) ;Y )/ exp 8 > < > : ~ S + t j 1 t ; ( ~ Y + t ) 1 + 1 9 > = > ; ~ S + t j 1 t ; ( ~ Y + t ) 1 (2.34) 39 The support of the above distribution is given by ~ S + t > 0; ~ Y + t 2 (l ; ; 1 2 ) for the posi- tive part and ~ S t < 0; ~ Y t 2 ( 1 2 ;l ; ) for the negative part. Following Buckle (1995) ([7]), for2 (1; 2) and2 (1; 1) we havet ; ( + 1 2 ) = + 1 andt ; (l ; ) = 0 hence p( + 1 2 j;X ( ~ Y + t ) ;Y ) =p(l ; j;X ( ~ Y + t ) ;Y ) = 0: Moreover, since t ; () is monotonic, it does not contribute to the maximum of (2.34). Hencep( ~ Y + t j;X ( ~ Y + t ) ;Y ) is a uni- modal distribution with mode at ~ Y + t =t 1 ; ~ S + t j 1 ! and the value at the maximum is p(t 1 ; ~ S + t j 1 ! j;X ( ~ Y + t ) ;Y ) = 1. This makes the rejection sampling algorithm suitable method to draw from the distribution in (2.34), for details on the rejection sampling please see Devroye 1986 ([17]) and Ripley (1987) ([43]). In this work, how- ever, we use adaptive rejection sampling algorithm (ARS) of Gilks and Wild (1992) ([24]), which utilizes the rejected draws 3 . This makes the sampling algorithm much more ecient especially in case of the distributions highly spiked around the mode. Since v + t = t ; ( ~ Y + t ) is a 1-1 correspondence, we obtain the draws on v + t and v t . In the derivation of other conditional distributions we condition onfv + t g T t=2 and fv t g T t=2 , which solves several multimodality problems as described below 4 . In the following we redene the vector of state variables by replacing ~ Y + t byv + t and ~ Y t by 3 We do not construct the envelope function as in Gilks in Wild (1992) ([24]) but follow closely Buckle (1995) ([7]). 4 The draw of v + t and v t for t = 2;::T contains information on the conditioning and parameters and hence it changes the property of updating procedure of parameter and if we condition on v + t and v t and not on ~ Y + t and ~ Y t . 40 v t for all t = 2;:::;T : X = (fh t g T t=1 ;f ~ U t g T t=2 ;f ~ S + t g T t=2 ;f ~ S t g T t=2 ;fv + t g T t=2 ;fv t g T t=2 ): (2.35) It is important to address the problem of 2 (1; 1) parameter. As noted above, we want to model the negative skewness of the asset returns, hence we have to consider > 0, moreover the degree of skewness of asset returns and work of Carr and Wu (2003) ([10]) suggests the choice of maximum negative skewness = 1 to be the reasonable one. However, in this setting we have to restrict to values 2 (0; 1), since according to proposition (1) we cannot guarantee unimodality of the distribution in eq. (2.34) for positive part ~ Y + t . The choice of bounds2 [0:1; 0:9] solves this problem and avoids over ow computer problems (for "too close" to 1 and "too close" to 0) in draws from the respective distributions. In the following we assume uniform prior distribution on 2 [0:1; 0:9]. The next problem involves choice of bounds for parameter 2 (1; 2). This is a delicate matter, since as & 1 the power coecient 1 in eq. (2.24) approaches innity. Moreover, as % 2 we approach normal distribution and loose identi- cation. Taking above into account we assume the uniform prior distribution on 2 [1:1; 1:99] to avoid the over ow computation problems, which was found to be not a restrictive assumption, since bounds are barely hit by the sampler. updating jump size variables ~ S + t and ~ S t for t = 2;:::T 41 By application of the Bayes rule the conditional posterior distribution for ~ S + t and ~ S t is given by: p( ~ S + t j ~ U t ;;X ( ~ S + t ) ;Y )/p(Y t jY t1 ;;X;Y (Yt;Y t1 ) )p( ~ S + t j ~ U t ;v + t ;;; j ) / exp ( 1 2 (Y t Y t1 (exph t1 ) 1 2 " (3) t S t ) 2 (exph t1 )(1 2 ) ) (2.36) exp 8 > < > : ~ S + t j 1 v + t 1 9 > = > ; ~ S + t 1 1 : By denition S t = ~ U t ~ S + t + (1 ~ U t ) ~ S t hence in the above density function the rst exponential part is a function of ~ S + t or ~ S t only if ~ U t = 1 or respectively ~ U t = 0. This property is intuitive, since there is no information contained in the sample about positive (negative) jump if there is a negative (positive) jump in the returns at time t. Please note, that the support of the density is ~ S + t > 0 for the positive part and ~ S t < 0 for the negative part. Let dene the following 1-1 correspondence: x + t = ~ S + t j 1 v + t : (2.37) Using the change of the variable formula, the density of x + t and x t is on the whole unimodal and log-concave moreover, this property is not aected if the sample con- tains information about the jump or not. Hence, we apply ARS algorithm by Gilks and Wild (1992) ([24]). We signicantly improve the ARS algorithm by supplying the unique maximum of the density forx + t andx t : In case, when the data contains no information about jump, it can be computed analytically by simple dierentia- 42 tion, in the other case we found the Newton's method to be ecient in computation. After draw of x + t we obtain the draw of ~ S + t by inverting the function in eq. (2.37). updating index of stability parameter 2 [1:1; 1:99] As noted in Buckle (1995) ([7]) updating the index of stability is the most dicult part in the Bayesian inference of stable jumps. We present his approach modied to accommodate for mixture of truncated stable distributions. He solved the problem of multimodality of complete conditional distribution by described above change of variable from the auxiliary variables ~ Y + t tov + t by the transformationv + t =t ; ( ~ Y + t ). If we condition not on n ~ Y + t o T t=2 and n ~ Y t o T t=2 but instead on v + t T t=2 and v t T t=2 the complete conditional distribution of is found to be on the whole unimodal and using Bayes rule is given by: p(j () ;X;Y )/p(f ~ U t g T t=2 j;)p(jf ~ S + t ;v + t g T t=2 ;f ~ S t ;v t g T t=2 ;; j )/ (2.38) p(f ~ U t g T t=2 j;)p(f ~ S + t ;v + t g T t=2 j;;; j )p(f ~ S t ;v t g T t=2 j;;; j )p(); wherep() is the prior distribution on; taken to be uniform on (1:1; 1:99);p() U(1:1; 1:99); independent of other parameters priors, and: p(f ~ U t g T t=2 j;) = (p ; ) P T t=2 ~ Ut (1p ; ) P T t=2 (1 ~ Ut) ; 43 p(f ~ S + t ;v + t g T t=2 j;;; j ) = 1 p ; T1 1 T1 exp 8 < : T X t=2 ~ S + t j 1 v + t 1 9 = ; T Y t=2 0 @ ~ S + t j 1 v + t 1 @t ; ( ~ Y + t ) @ ~ Y + t 1 t ; ( ~ Y + t )=v + t 1 A ; p(f ~ S t ;v t g T t=2 j;;; j ) = 1 1p ; T1 1 T1 exp 8 < : T X t=2 ~ S t j 1 v t 1 9 = ; T Y t=2 0 @ ~ S t j 1 v t 1 @t ; ( ~ Y t ) @ ~ Y t 1 t ; ( ~ Y t )=v t 1 A : Note, that the state vector X is already redened in eq. (2.35) and contains in- formation on v + t T t=2 and v t T t=2 . In order to compute the value of the above conditional distribution in eq. (2.38) we need to inverse function t ; for given v + t and v t and nd the respective values of ~ Y + t and ~ Y t . This can be done eciently using Newton's method as suggested in Buckle (1995) ([7]). Since the complete con- ditional distribution of is of unknown form, we rely on the Metropolis-Hastings (MH) step to sample from this distribution. The random walk MH step with normal proposal distribution was found to be ecient. updating skewness parameter 2 [0:1; 0:9] 44 We take independent uniform prior on;p()U(0:1; 0:9). Updating skewness parameter is similar to the update of : p(j () ;X;Y )/p(f ~ U t g T t=2 j;)p(jf ~ S + t ;v + t g T t=2 ;f ~ S t ;v t g T t=2 ;; j )/ (2.39) p(f ~ U t g T t=2 j;)p(f ~ S + t ;v + t g T t=2 j;;; j )p(f ~ S t ;v t g T t=2 j;;; j )p(); where all components are derived in the updating procedure for but should be treated as the functions of . By conditioning on v + t T t=2 and v t T t=2 we achieve unimodality of complete conditional posterior of. The random walk MH step with normal proposal distribution is also found to be ecient in this case. updating scale parameter j > 0 Let assume the following inverse gamma prior on j ,p( j )IG(dd;DD); where dd = 3 andDD = 1 20 ; which is relatively uninformative prior. We have the following complete conditional posterior for j : p( j j ( j ) ;X;Y )/p(f ~ S + t ;v + t g T t=2 j;;; j )p(f ~ S t ;v t g T t=2 j;;; j )p( j )/ (2.40) exp 8 < : 1 j 1 2 4 T X t=2 0 @ ~ S t 1 v t 1 + ~ S + t 1 v + t 1 1 A 3 5 9 = ; 0 @ 1 1 j 1 A 2(T1) p( j ): Let dene the proposal density in the MH step for updating j to be: q( j )/ exp 8 < : 1 j 1 2 4 T X t=2 0 @ ~ S t 1 v t 1 + ~ S + t 1 v + t 1 1 A 3 5 9 = ; 0 @ 1 1 j 1 A 2(T1) : 45 We can directly draw from this distribution by the change of variable. Let ^ j = 1 j , then by change of variable formula, we have the following density for ^ j : ^ q(^ j )/ exp 8 < : 1 ^ j 2 4 T X t=2 0 @ ~ S t 1 v t 1 + ~ S + t 1 v + t 1 1 A 3 5 9 = ; 1 ^ j 2(T1) 1 ^ j 1 : We recognize in this density the inverse gamma density, hence ^ j IG((2(T 1) + 1 1);DD ); where DD = T X t=2 0 @ ~ S t 1 v t 1 + ~ S + t 1 v + t 1 1 A : After drawing from this distribution, we have to solve for j = ^ 1 j . The MH algorithm with the above choice of proposal density is very ecient with acceptance rate of above 98%. This also re ects that our prior is uninformative, since for completely at prior on j the acceptance rate should be 100%. updating auxiliary mixing variables ~ U t for t = 2;:::T Sincep( ~ U t j;) is a Bernoulli distribution with probability of successp ; , hence the complete conditional posterior is also Bernoulli and is given by: p( ~ U t = 1j;X ( ~ Ut) ;Y )/p(Y t j;X ~ Ut ; ~ U t = 1;Y t1 )p( ~ U t = 1j)/ (2.41) exp ( 1 2 (Y t Y t1 ~ S + t (exph t1 ) 1 2 " (3) t ) 2 exph t1 (1 2 ) ) p ; ; and analogously: p( ~ U t = 0j;X ( ~ Ut) ;Y )/p(Y t j;X ( ~ Ut) ; ~ U t = 0;Y t1 )p( ~ U t = 0j)/ (2.42) exp ( 1 2 (Y t Y t1 ~ S t (exph t1 ) 1 2 " (3) t ) 2 exph t1 (1 2 ) ) (1p ; ): 46 We can directly calculate the above probabilities and normalize their sum to unity. The draw from this distribution is then straightforward. 2. Model SJSV (2.20) It is important to note, that update procedure of SJSV model jump parameters and jump state variables is similar to that of models SVSJ and SVSJLE. In the above it is enough to substitute 2 for exph t for allt and = 0 in order to achieve constant volatility from diusion and x value of j = 1, since this parameter no longer is estimated and is replaced with stochastic volatility. Moreover, in conditional posteriors in eq. (2.36), (2.41) and (2.42) we have to premultiply all jump size values S t ; ~ S + t and ~ S t respectively by exp( t1 ) value for all t = 2;:::;T . updating stochastic volatility states t from Levy stable jump part for t = 1;:::T . For t = 1 we have: p( 1 j;X ( 1 ) ;Y )/p(Y 2 jY (Y 2 ) ;;X)p( 1 j 2 ; ; ; ); where the second probability is given by the symmetry formula for autoregressive models (AR): p( 1 j 2 ; ; ; )N( + (1 ) 2 ; 2 ); 47 and the rst is given by: p(Y 2 jY (Y 2 ) ;;X)/ exp( 1 2 (Y 2 Y 1 exp( 1 )S 2 ) 2 2 ): In the similar manner we can derive the complete conditional posteriors for 1<t< T : p( t j;X (t) ;Y )/p(Y t+1 jY (Y t+1 ) ;;X)p( t j t+1 ; t1 ; ; ; ); where p( t j t+1 ; t1 ; ; ; )N( (1 )( t+1 + t1 ) + ( ) 2 (1 ) 2 + 1 ; 2 (1 ) 2 + 1 ); and p(Y t+1 jY (Y t+1 ) ;;X)/ exp( 1 2 (Y t+1 Y t exp( t )S t+1 ) 2 2 ): After careful study of conditional posteriors for t for t = 1; 2;:::; (T 1) we found, that they are either bimodal or unimodal and in general in the unimodal case they are not log-concave. However, we can numerically compute all local max- imas using Newton's method without much computational burden and then apply rejection algorithm utilizing the rejected points to update the envelope log-density function. Finally, the complete conditional posterior for T is N( + (1 ) T1 ; 2 ). updating > 0 48 In model SJSV we incorporate stochastic volatility from jumps only and take the constant volatility > 0 from diusion. We take almost uninformative inverse gamma prior on 2 given by 2 IG(d;D), withd = 3 andD = 1 20 . By the standard OLS theory the complete conditional posterior of 2 is also inverse gamma. 2.5 Auxiliary Particle Filter In the following, we restrict our attention to models with stochastic volatility from diusion 5 . In order to perform value at risk (VaR) analysis we x vector of parame- ters for each model at the respective MCMC estimate and calculate the following probabilities: Pr(R t+L <r t+L j(m);R t ;m); (2.43) L2f1; 5g; t2N L ; R t =fR s g s=1;:::;t m2fSVPJ;SVPJLE;SVVG;SVVGLE;SVSJ;SVSJLEg where (m) denotes parameter vector for model m specication, R t+L = Y t+L Y t+L1 is daily log-return on the asset at time t +L with its law determined by the model specication and r t+L is observed value of this log-return at time t +L. In the above N L denotes the subset of natural numbers less than TL and divisible 5 We do not develope auxiliary particle lter procedure for models with stochastic volatility from jumps, since in section (3.3) we found them to perform the worst in terms of asset returns representation . 49 byL. This guarantees, that we analyze only "non-overlapping" periods, but reduces the number of available data to roughly one-fth in case of L = 5. In this work we consider one-day ahead and ve-day ahead time horizons for VaR analysis, which makes possible to assess model ability to predict one-day ahead and one-business week ahead daily log-returns. Note, that (2.43) can be calculated not only for the in-sample period but also for out-of sample period, whenever we have data available. We can study VaR performance of the model by comparison of given signicance level and unconditional covering frequency of each model implied by covering probabilities in (2.43). Moreover, if model is correctly specied such binary variables, indicating if data points are contained in the VaR interval, should be independently distributed and hence there should be no "clustering" in time of their respective realizations. We can estimate values in eq. (2.43) by: z (L) t = c Pr(R t+L <r t+L j;R t ;m) = (2.44) 1 S S X s=1 Pr(R t+L <r t+L j;R t ;m;h (s) t+L1 ;h (s) t+L ;J (s) t+L ); where (h (s) t+L1 ;h (s) t+L ;J (s) t+L )iid p(h t+L1 ;h t+L ;J t+L j;R t ;m). The draws for L = 1 can be performed by utilizing following condition: p(h t ;h t+1 ;J t+1 j;R t ;m) =p(h t+1 jh t ;;m)p(h t j;R t ;m)p(J t+1 j;m): (2.45) 50 Things complicate forL = 5, where we have to draw from the following joint distri- bution: p(h t+4 ;h t+5 ;J t+5 ;h t+3; h t+2 ;h t+1 ;h t j;R t ;m) (2.46) which can be further decomposed to: p(h t+5 jh t+4 ;;m)p(h t+4 jh t+3 ;;m)p(h t j;R t ;m)p(J t+5 j;m): By discarding draws for variables, that do not directly enter in equation (2.44), we have draws from the desired p(h t+4 ;h t+5 ;J t+5 j;R t ;m) distribution. It is important to note, that above equations hold for all models with stochastic volatility from diusion with jump sizes: J t (m) = 8 > > > > > > < > > > > > > : q t t for models with Poisson jumps J VG t for models with VG jumps S t (;; 0; j 1= ) for models with stable jumps If given model is correctly specied, then z (L) t (for t2 N L ) should be iid U(0; 1) distributed. By transformation using inverse cumulative distribution function of standard normal distribution, F 1 , we dene: b z (L) t =F 1 (z (L) t ): (2.47) The distribution ofz (L) t implies, that transformed variablesb z (L) t should beiidN(0; 1) distributed. This fact is later used for model evaluation in view of its predictive performance. 51 In order to sample from distributions in eq. (2.45) and (2.46) we have to sample from ltering densityp(h t j;R t ;m) and then, conditional on this draw, sample from all predicting densities p(h t+l jh t+l1 ;;m) and p(J t+L j;m). Sampling from these densities is straightforward and the most dicult problem involves approximation of ltering density p(h t j;R t ;m) by auxiliary particle lter developed in Pitt and Shephard (1999) ([42]). Chib, Nardari and Shephard (2002) ([15]) extend basic auxiliary particle lter of Pitt and Shephard (1999) ([42]) for Poisson type jumps in returns but do not include leverage eect. Durham (2004) ([18]) extends the basic particle lter for models with leverage eect but does not include jumps in returns, moreover, he works with particle lter and does not apply auxiliary particle lter involving index parameters draws 6 . In this work, we ll this gap in the literature by presenting auxiliary particle lter for models with both: jumps in returns and leverage eect. Let rst notice, that: p(h t+1 ;J t+1 (m)jm;;R t+1 ) = Z p(h t+1 ;h t ;J t+1 jR t+1 ;) p(h t j;R t ) dp(h t j;R t ) (2.48) and p(h t+1 ;h t ;J t+1 jR t+1 ;) = p(h t+1 ;h t ;J t+1 ;R t+1 jR t ;) p(R t+1 jR t ) = (2.49) p(R t+1 jh t+1 ;h t ;J t+1 ;R t ;)p(h t+1 jh t ;)p(J t+1 j)p(h t j;R t ) p(R t+1 jR t ) 6 For discussion on improvements of auxiliary particle lter over particle lter please refer to Pitt and Shephard (1999) ([42]). 52 where we omitted model specication m for notational simplicity. Substituting (2.49) into (2.48) we have: p(h t+1 ;J t+1 j;R t+1 )/ (2.50) Z p(R t+1 jh t+1 ;h t ;J t+1 ;R t ;)p(h t+1 jh t ;)p(J t+1 j)dp(h t j;R t ) Auxiliary particle lter is a recursive algorithm to approximate ltering densi- ties p(h t j;R t ) for t = 0;:::;T by nite number M of "particles" for each t, that dene discrete probability distribution lter b p(h t j;R t ;m). We denote particles for lter at time t as h (k) t , where k = 1; 2;:::;M. Given M particles dening dis- crete probability distribution lterb p(h t j;R t ;m) at timet, we obtain approximation b p(h t+1 j;R t+1 ;m) for t + 1 dened by its respective M particles using relation in (2.50): 1. draw R M indexes k 1 ;k 2 ;:::;k R from discrete distribution g(kjR t+1 ) with support of k = 1; 2;:::;M; with weights proportional to: g(kjR t+1 )/p(r t+1 j k t+1 ;h k t ; b J t+1 = 0;R t ;) (2.51) Note, that k represents index on mixture in (2.50) as in Pitt and Shep- hard (1999) ([42]). We x b J t+1 = 0 and k t+1 to be the mean of prediction density p(h t+1 jh k t ;), which is straightforward to calculate for all model spec- ications. Denote index draws as b k 1 ; b k 2 ;:::; b k R and record respective values of b k 1 t+1 ; b k 2 t+1 ;:::; b k R t+1 and h b k 1 t ;h b k 2 t ;:::;h b k R t . 53 2. Draw proposal particles h (s) t+1 ; s = 1; 2;:::;R given mixture index and particle from preceding lter using respective prediction density: h (s) t+1 p(h t+1 jh b ks t ;) 3. Draw jumps J (s) t+1 ; s = 1; 2;:::;R from p(J t+1 j;m), where we include index model m for three dierent jump specications. 4. Reweight draws (h (s) t+1 ;J (s) t+1 );s = 1; 2;:::;R by drawingM times (with replace- ment) from discrete probability distribution with weights proportional to: $ s = p(r t+1 jh (s) t+1 ;h b ks t ;J (s) t+1 ;R t ;) p(r t+1 j b ks t+1 ;h b ks t ; b J t+1 = 0;R t ;) (2.52) for s = 1; 2;:::;R. We nally get M draws dening discrete lter distribution b p(h t+1 j;R t+1 ;m) by discarding draws on J t+1 . 5. Go to point 1. and increment t. In our applications we take R = 30; 000 and M = 20; 000. For L = 1 and L = 5 we choose respectively S = M and S = 10M in eq. (2.44). In case of L = 1 we do not draw from discrete auxiliary particle lter distributionb p(h t j;R t ;m) but directly utilize all derived particles. The above choice of S;R and M constants was sucient to induce low variability of statistics calculated using derived particles among dierent starting seeds of random number generator. 54 Please note, that above algorithm works for all model specications with stochas- tic volatility from diusion. If model does not include leverage eect all of the above applies and can be further simplied. 55 Chapter 3 Estimation Results 3.1 The Data In this work the data on the S&P 500 index extends from 01/02/1981 to 08/31/2007 and comprises 6730 daily observations available for download atfinance:yahoo:com The S&P 500 index levels are reported at the closing times in each business day. The data is divided into two subsamples for the in-sample period covering 01/02/1981- 02/02/2007 and out-of sample period covering 02/02/2007-08/31/2007. The eight model specications are estimated using the in-sample period data, that gives a total of T = 6584 daily observations on the S&P 500 index. This allows for modeling market crash in 1987 and the \dot.com" corrections from 1999 to 2001. The out-of sample period of 146 observations is left for the value at risk (VaR) analysis and used as one of the models selection criteria. This allows for checking the model out-of sample performance in the summer 2007 corrections. The modeling of S&P 500 index returns was already extensively addressed in the literature, since its deviation from assumptions behind the standard Black- 56 Scholes (1973) ([6]) model are undisputable. In table (1) and gure (1) we present respectively the descriptive statistics of daily log-returns on the S&P 500 index and graphs of S&P 500 index log-level and S&P 500 index log-returns. Since the log-returns in the Black-Scholes model are normally distributed, there is no place for modeling of excess kurtosis and skewness. The data on the S&P log-returns, however, indicate that there exists signicant negative skewness of1:7927 and kurtosis of 44:2079 1 . The out-of sample period characteristics dier in all respects from the in-sample period. The ability of each model to predict these behavior of S&P 500 index returns is therefore of much importance. In gure (2) we also report the skewness and kurtosis as a term structure of S&P 500 returns. The speed of convergence towards normal distribution is also one of the important characteristics of the data. Table 1. S&P500 daily log-returns descriptive statistics. 1 This values are reported for the in-sample period. 57 Figure 1. S&P 500 index log-level and log-returns. S&P 500 index log-level and log-returns, daily data 01/02/1981-08/31/2007. 58 Figure 2. Skewness and kurtosis of S&P500 returns for term structure of returns. Data on returns 01/02/1981-02/02/2007. We consider 1-day up to yearly horizon returns. 59 3.2 Identication problem By introduction of innite activity jumps in returns we would obviously achieve more degree of freedom in representing the desired behavior of the returns or any other series under consideration. This generalization, however, leaves the question if discretely sampled data contains enough information about this type of jumps and if we are able to disentangle them from the diusion part, and if the answer is positive, what sample size is enough to guarantee the required degree of estimation precision. The above is found not to be a problem in the estimation under risk neutral measure, since the options data contains much more information on the tails of the underlying asset returns process. The estimation results of the time changed returns process with stable jumps or VG jumps under risk-neutral measure can be found in Huang and Wu (2004) ([27]), where authors apply the fast Fourier transform (FFT) method to nd the model implied option prices to be later compared with the observed option prices. In this work we are interested in the estimation under statistical measure, which can only utilize information on the observed asset returns. The recent work by A t-Sahalia (2003) ([1]) provides positive theoretical answer for the simple returns model with Cauchy jumps (stable jumps with = 1 and = 0) and with constant volatility from diusion. We need to show that this property also holds for 2 (1:1; 1:99), and for negative skewness 2 (0:1; 0:9) 60 and in presence of stochastic volatility. It is important to note, that for = 2 the importance of skewness cancels out and we are left with normal distribution. Hence, as approaches two we need more information to distinguish between stable jumps and diusion in the model without stochastic volatility, since most impor- tant information is contained in the very tails of the returns distribution. In our specication, we introduce stochastic volatility in the form of log-volatility, which was found to improve the model ability to disentangle time changed diusion (dif- fusion with stochastic volatility) from stable jumps in models SVSJ and SVSJLE even in case of moderate sample sizes. Moreover, we also veried identication in model SJSV in presence of constant volatility from diusion and time changed stable jumps (stochastic volatility from jumps). In case of models with latent VG jumps with stochastic volatility from diusion, Li, Wells and Yu (2006) ([36]) empirically verify identication. In our work, however, we also have to deal with stochastic volatility from VG jumps and this problem under statistical measure, to our best knowledge, is not yet addressed in the literature. Since empirical verication of identication involves simulation of many sample paths from data generating process (DGP) for given model specication and param- eter values, it requires signicant computational resources. In our work we decided to choose specic parameter values for the DGP for each model. At rst we run MCMC sampler and found parameter estimates for the S&P 500 data for all models 2 2 Please refer to section (3.3) for details. 61 provided in table (3). Then, we reported posterior means as the parameter estimates and used them for the DGP specication. Hence, in our reverse-engineering experi- ment we check for identication only for the parameters, that are found to be MCMC estimates for S&P 500 data. Moreover, we generate the same number of observations T = 6584 in sample paths from each DGP as the sample size of S&P 500 data used for estimation. This would allow to answer the question if such sample size is enough to estimate the parameters with high degree of precision. Optimally we would like to generate large number of sample paths from such specied DGPs, however, in our work we restricted to only one sample path from each model taken with the starting seed of the random number generator inMatlab 7:0:4. We leave bias and root mean square error (RMSE) calculation based on large number of generated sample paths for future research. Moreover, we do not perform identication analysis for models SVVG and SVSJ, since they can be seen as special cases of models SVVGLE and SVSJLE respectively with constraint = 0. Let consider model specications with VG jumps in returns SVVGLE, VGSV respectively in eq. (2.16) and (2.17). In each case we treat rst 1; 000; 000 draws from the MCMC sampler as the burn-in draws and treat the next 1; 000; 000 as the draws from the posterior distribution. The starting values of the mean/drift parameters ; h ; ; j are chosen to be zero, the scale parameters h ; ; j ; to 62 be one, the correlation parameter to be zero and h ; and j to be one 3 . The starting values for the latent variables involve only choice of jump sizes J VG(1) t = 0 and (1) t = 0; h (1) t = 0 for all t. Since we draw the latent variables G t at the very beginning of the MCMC algorithm, we do not need to specify these values. We calculate the mean and the standard deviation of the posterior distributions for each of the parameters and treat the former as the parameter estimate. The parameter estimates were found not to be aected by dierent choices of starting values for the MCMC algorithm. Descriptive statistics of posterior distributions for simulated data are presented in table (2). All parameters are estimated with reasonably high degree of precision even with moderate sample size of T = 6584 observations. In the next step we perform similar analysis for models with Levy stable jumps in returns SVSJLE and SJSV respectively in eq. (2.19) and (2.20). We treat rst 100; 000 draws from MCMC chains for each model as the burn-in period and the next 100; 000 as the draws from stationary distribution. The starting values for mean/drift parameters; h ; are zero, the scale parameters h ; ; j ; are one, the correlation parameter is zero, h ; are one and jump specic parameters (1) = 1:5 and (1) = 0:5 . The choice of starting values for the latent variables involves choice of ~ U (1) t = 1 (only positive jumps), h (1) t = 0, (1) t = 0 and nally ~ S +(1) t = 0:01 and ~ S (1) t =0:01 for all t. We present the parameter estimates for 3 Note, that only some of those pramaters refer to each model, if some parameter is not in the DGP specication, it should be disregarded. 63 each of the models in table (2), which are found not to be aected by dierent choices of starting values of MCMC algorithms. All parameters are identied in models with Levy stable jumps and estimated with reasonably high degree of precision with low posterior standard deviations. Summing up, the discrete data contains enough information to distinguish in- nite activity jumps with either nite variation (VG jumps) or innite variation (Levystable jumps) from diusion part of the model in the presence of stochas- tic log-volatility introduced both via diusion or jump part. The sample size of 6584 observations is enough to achieve high degree of estimation precision. 64 Table 2. Descriptive statistics of posterior distributions for models with innite activity jumps for simulated data. For each model specication we simulate sample of T = 6584 observations from respective true DGP with parameters xed at the MCMC estimates for S&P500 data. 65 3.3 Estimation Results for S&P 500 Data As described in section (3.1) we have T = 6; 584 daily observations on the S&P 500 index level used for parameter estimation. In table (3) we present descriptive statistics of posterior distributions for all model specications. Table 3. Descriptive statistics of posterior distributions for S&P500 return models. Mean and standard error (in parentheses) of respective posteriors, S&P500 index data on T = 6584 observations, 01/02/1981-02/02/2007. 66 The MCMC draws along with smoothed log-volatilites and jump sizes are illus- trated in gures (3) - (10) for all eight models depending on the source of stochastic volatility and jump structure. The starting values of the parameters, burn-in pe- riods and number of draws from the stationary distributions for models with VG jumps and Levy -stable jumps are the same as for simulated data and given in section (3.2). For models with compound Poisson jumps in returns the drift and log-volatility related parameters and latent variables are given the same starting values as for models with innite activity jumps. Poisson jump specic parameters are given starting values (1) j = 0, (1) j = 1 and (1) j = 0:5 and for latent variables we assume no jumpsq (1) t = (1) t = 0 for allt. The burn-in period is taken to be 100; 000 draws and the next 100; 000 are treated as draws from stationary distribution. In the following we also implement the model selection criteria developed in Jones (2003) ([33]). From assumption, in model specications (2.13)-(2.20) " (1) t and " (2) t are iid N(0; 1). In the following let call " (1) t as residuals from returns equation, " (2) t as residuals from log-volatiltiy equation in models without leverage eect and nally " (1) t + p 1 2 " (2) t as residuals from log-volatiltiy equation in models with leverage eect. We may view those residuals as latent variables. Hence, we can construct posterior distributions for functions of these latent variables by evaluating these functions at each step of the MCMC algorithm. Since model residuals are assumed to beiidN(0; 1) we calculate mean, standard deviation, skewness, kurtosis 67 Figure 3. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (1) specication with stochastic volatility and no leverage eect, compound Poisson jumps. 68 Figure 4. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (2) specication with stochastic volatility and leverage eect, compound Poisson jumps. 69 Figure 5. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (3) specication with stochastic volatility from diusion, variance gamma jumps and no leverage eect. 70 Figure 6. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (4) specication with stochastic volatility from diusion, leverage eect and variance gamma jumps. 71 Figure 7. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (5) with stochastic volatility from variance gamma jumps. 72 Figure 8. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (6) specication with stochastic volatility from diusion, stable jumps and no leverage eect. 73 Figure 9. MCMC draws for S&P 500 data (01/02/1981-02/02/2007), model (7) specication with stochastic volatility from diusion, leverage eect and stable jumps. 74 Figure 10. MCMC draws for S&P500 data (01/02/1981-02/02/2007) model (8) with stochastic volatility from stable jumps. 75 and rst-order autocorrelation. Then we calculate the median and 95% condence intervals for those statistics reported in table (4). A correct model specication implies, that mean is zero, standard deviation one, skewness zero, kurtosis three and autocorrelation zero. The important feature of our analysis is, that estimation of all eight model spec- ications allows for drawing conclusions about the marginal importance of leverage eect, dierent jump structures and source of stochastic volatility. 1. Leverage Eect According to importance of leverage eect specication, we closely follow Jacquier, Polson and Rossi (2004) ([29]). Since we enjoy the same log-volatility specication, their results are applicable in our case. In their work, the inclusion of leverage eect is found to correct for possible misspecications resulting in biased estimates of volatility states and log-volatility process parameters. Models with Poisson jumps, SVPJ and SVPJLE Inclusion of leverage eect corrected log-volatility parameters estimates resulting in higher estimate of Poisson intensity j in model with leverage eect (change from 0:0016 to 0:002) and closer to zero j jump mean (change from0:0537 to0:049). The parameter controlling speed of mean reversion of log-volatility is corrected away from non-stationarity bound for model with leverage eect. In both models the least 76 Table 4. Posterior medians and 95% condence intervals. Posterior medians and 95% condence intervals (in parentheses) of mean, std. deviation, skewness, kurtosis and daily autocorrelation of model residuals. 77 precise estimates involve mean jump parameters j . This fact can be explained by small number of realized jumps in returns implied by low values of Poisson intensity j . Introduction of leverage eect improves model ability to represent both skewness and kurtosis of returns distribution as presented in table (4). Both models are estimated with high degree of precision. The leverage eect in model SVPJLE is found to be signicant with estimate of =0:5767. Models with VG jumps, SVVG and SVVGLE In case of models with VG jumps in returns the leverage eect has a dramatic eect on the parameter estimates and precision of estimation. The serious misspec- ication of the model without leverage eect is documented by hitting the non- stationarity bounds for h parameter in the MCMC algorithm with draws in gure (5) and therefore imprecise estimation of all log-volatility parameters. Hence, inclu- sion of leverage eect in models with VG jumps in returns is of special importance. The leverage eect is estimated at0:7928 with high degree of precision and esti- mation of models with zero constraint on is a serious source of misspecication. In both models the jump specic parameters are estimated with high degree of precision with an exception of j parameter controlling skewness of returns. Models with stable jumps, SVSJ and SVSJLE Inclusion of leverage eect for models with stable jumps in returns has standard eect of correction of log-volatility equation estimates. However, this impact is 78 not as dramatic as in case of models with VG jumps in terms of both estimates and estimation precision. We also witness correction of h parameter estimate away from non-stationarity bound for model with leverage eect. The most important change involves skewness parameter estimates from 0:234 for model without leverage eect to 0:4664 for model with leverage eect. This can be explained by the ability of representing skewness of returns by both leverage eect and skewness parameter value . The leverage eect coecient is estimated at level0:7605 with high degree of precision. We also observe slight improvement in tting returns equation documented in table (4) for model with leverage eect. 2. Jump structures for models with stochastic volatility from diusion In this analysis we compare models with leverage eect only, since in models with stochastic volatility from diusion the leverage eect is found to have an im- portant role in explaining S&P500 returns behavior. We consider models SVPJLE, SVVGLE and SVSJLE in our comparison. The problem in estimation precision involves parameter j for VG model, j for Poisson jumps model. Since all other parameters are estimated with high degree of precision reported in table (3), we focus our attention on residuals characteristics from the respective returns and log- volatility equations in table (4). Since Poisson type jumps occur very rarely with estimate of Poisson intensity of j = 0:0020, which gives about one jump per two years, they are unable to t the characteristics of the return process. The "small" 79 and frequent jumps are simply not represented by model with Poisson type jumps. As we can expect the skewness, aected by occurrence of large jumps in the very tails of the distributions, is much better represented, than kurtosis of returns as presented in table (4). On the other side, models with innite activity jumps in returns, that is the models with innite number of "small" jumps in nite time interval, have almost perfect t with respect to returns and log-volatility equations. Since the dierences between skewness and kurtosis representation are negligible, with slight advantage of Levy -stable model, we rely on the VaR analysis as model selection criteria presented in section (3.4). 3. The source of stochastic volatility Models VGSV and SJSV with respectively VG jumps and Levystable jumps in returns and stochastic volatility from jumps are estimated with high degree of pre- cision presented in table (3). The lack of stochastic volatility from diusion dwarves both models ability to represent returns characteristics especially with respect to kurtosis of approximately 3:6 value for both models reported in table (4). Their performance is even worse, than that of the models with Poisson type jumps only. Most importantly, we nd, that diusion is still an important part of all models even with innite activity jumps and innite total variation jumps, although their 80 characteristics resemble much closer behavior of the Brownian motion, than the rare jump Poisson type process. 3.4 VaR analysis In this section we apply auxiliary particle lter described in section (2.5) to evaluate one-day ahead and one-week ahead prediction and VaR performance of models with stochastic volatility from diusion. We do not consider models with stochastic volatility from jumps, since in section (3.3) we found them to perform the worst in terms of asset returns representation. In tables (5) and (6) we present descriptive statistics offb z (L) t g distribution for L = 1 and L = 5 respectively, withb z (L) t dened in eq. (2.47) calculated for dierent model specications with their respective parameters xed at the MCMC estimates as in table (3) for S&P500 data on period 01/02/1981-02/02/2007. Data on S&P500 daily returns, used forb z L t calculation, are derived from S&P500 data for all available observations from 01/02/1981-08/31/2007. A correct specication implies, that fb z L t g is iid N(0; 1) distributed, hence the mean is zero, standard deviation one, skewness zero, kurtosis three and 1st. order autocorrelation zero. For one-day ahead analysis (L = 1) in table (5) models with leverage eect pro- duce slightly decreased prediction performance comparing to models with the same jump structure but no leverage eect across all specications with innite activity 81 Table 5. Descriptive statistics offb z (1) t g distribution. Table 6. Descriptive statistics offb z (5) t g distribution. 82 jumps in returns. It is important to note, however, that this decrease is not signi- cant in case of Levy-stable jump models, and in the view of results from previous section, we avoid serious misspecication problems by incorporating leverage eect. For models with Poisson jumps inclusion of leverage eect does not have such char- acteristics for one-day ahead prediction performance. Most importantly, we nd that models with VG jumps in returns perform the worst in terms of one-day ahead prediction, while models with Poisson and Levy -stable jumps produce satisfac- tory results. Although models with VG jumps represent skewness and kurtosis of asset returns well as presented in section (3.3), one source of potential problems is re ected in lack of j estimation precision controlling skewness of returns and hence tails of daily returns distribution. It is important to note, that values b z (1) t depend on the left tail of the one-day ahead conditional returns distribution and its skew- ness also depends on j . On the other side, since models with Poisson jumps also suer from lack of precision estimation of j parameter also controlling skewness of returns distribution, its eect is not as strong, since jumps occur very rarely with small estimates of Poisson jump intensity j . For one-week ahead analysis (L = 5) in table (6) models with VG jumps improve signicantly with respect to prediction performance in view of skewness and kurtosis but also completely fail in terms of standard deviation, which departs signicantly from unity. In case of models with Poisson jumps and Levy-stable jumps increase 83 in time horizon from one-day to one-week decreased performance in terms of kurtosis but corrected for degree of skewness. Improvement in terms of serial correlation across all model specications is natural, since we take longer time horizon into consideration. In tables (7) and (8) we present one-day ahead (L = 1) and one-week ahead (L = 5) VaR performance of the models respectively. We compare signicance levels of 1%, 5% and 10% to respective empirical coverage probabilities for each model. We calculate empirical coverage probabilities by computing mean of binary variables representing visits of observed daily returns to the model implied VaR regions. This is done with use ofz (L) t variables in eq. (2.44), which are a by-product of computations described above. We separate analysis for the in-sample period 01/02/1981-02/02/2007 and out-of sample period 02/03/2007-08/31/2007 for one- day ahead analysis and consider the whole sample for one-week ahead analysis. The VaR analysis for one-week ahead time horizon cannot be separately performed for out-of sample period due to small number of observations. We nd, that inclusion of leverage has an eect of increase (in absolute value) of VaR estimates and this property does not depend neither on the jump structure nor on the time-horizon. Moreover, models with VG jumps in returns perform poorly in the very left tail of conditional returns distribution for all signicance levels and both time-horizons L = 1 and L = 5. They systematically overestimate VaR values (in 84 Table 7. One-day ahead value at risk analysis (VaR) performance of models with stochastic volatility from diusion. Table 8. Five-day ahead value at risk analysis (VaR) performance of models with stochastic volatility from diusion. 85 absolute value), which results in too small covering probabilities for observed data points. On the other side, models with Poisson and Levy -stable jumps perform well for both one-day and one-week ahead. Interestingly enough VG jump models outperform Poisson and Levy -stable jump models in one-day ahead analysis for out-of sample period. However, due to very small number of available observations, we cannot draw far reaching conclusions from this result and need to extend analysis when more data become available. Summing up, the only specications, that represent well one-day and one-week ahead conditional distribution of returns and their left tails are models with Poisson jumps and Levy stable jumps. However, since models without innite activity jumps and leverage eect do not represent data well and can be source of potential misspecication as described in section (3.3) we conclude, that model with Levy -stable jumps and leverage eect dominates other specications. 86 Chapter 4 Conclusions In this work we use data on daily S&P500 index returns and perform estimation under statistical measure. We consider model specications with dierent jumps in returns specications, dierent source of stochastic volatility and either existence or lack of leverage eect. The source of stochastic volatility is taken to be either from diusion or from the jump part in case of models with innite activity jumps. This enables us to analyze the marginal eects of jumps, source of stochastic volatility and leverage eect. We face the problem of parameters estimation in presence of both latent volatility and latent jump sizes. We perform estimation using Monte Carlo Markov chain methods and propose special algorithm for models with Levy-stable jumps, which solves the problem of MCMC state-space separability. The models with innite activity jumps in returns (both VG and stable jumps) are able to represent well excess kurtosis of the underlying data and tails of the returns distribution if diusion is the source of stochastic volatility. Models with stochastic volatility from pure jump part does not represent data well and this fea- 87 ture does not depend on the inclusion of the leverage eect from diusion, since models with stochastic volatility from diusion but no leverage eect also perform signicantly better. Moreover, the leverage eect is found to be an important com- ponent of all models. In case of models with innite activity jumps it has an eect of higher degree of precision of parameter estimates by reducing the misspecication of independent increments of returns and log-volatility. In order to evaluate one-day and one-week ahead prediction performance of the models we construct an auxiliary particle lter algorithm for models with stochastic volatility from diusion, leverage eect and jumps. We nd, that models with VG jumps in returns perform the worst in terms of one-day and one-week ahead pre- diction, while models with Poisson and Levy stable jumps produce satisfactory results. This inability of VG jump models to represent one-day and one-week ahead conditional returns distribution can be explained by the lack of estimation precision of parameter governing skewness of returns. Value at risk (VaR) analysis, which focuses on the left-tail characteristics of conditional returns distribution, also con- rms poor performance of models with VG jumps. Models with Poisson jumps and Levy stable jumps perform satisfactory. In general, VG jump models tend to overestimate VaR values (in absolute value). Finally, an inclusion of leverage eect from diusion is found to increase VaR estimates (in absolute value) across all jump specications. 88 The above results suggest, that models with stochastic volatility from diusion and Levy stable jumps in returns represent well distribution of returns and can be used for prediction and VaR analysis. A line for future research is to study the importance of jumps in volatility simi- larly to Eraker, Johannes and Polson (2003) ([20]). Their specications incorporate only Poisson type jumps and diusion as the only source of stochastic volatility, hence extensions involving innite activity jumps and jumps as source of stochastic volatility is of much importance. Moreover, in this work we analyze diusion as the only source of leverage eect and further research could involve model specications with dierent source of leverage eect in presence of innite activity jumps and stochastic volatility from jump part. This is an interesting problem, since innite activity jumps resemble diusion behavior much closer, than rare and "big" Poisson jumps. 89 References [1] A t-Sahalia, Y., 2003, Disentangling volatility from jumps, Working Paper [2] Andersen, T.G., L. Benzoni and J. Lund, 2001, An empirical investigation of continuous-time equity return models, Working Paper no. 8510 (National Bureau of Economic Research) [3] Applebaum, D., 2004, Levy processes and stochastic calculus (Cambridge Uni- versity Press). [4] Bates, D.S., 1995, Post-crash moneyness biases in S&P 500 futures options, working paper (Rodney, L. White, Wharton School, University of Pennsylvania, Philadelphia, PA.) [5] Bertoin, J., 1998, Levy processes (Cambridge University Press, Cambridge) [6] Black, F. and M. Scholes, The pricing of options and corporate liabilities, Jour- nal of Political Economy, 81, 637-659. [7] Buckle, D.J., 1995, Bayesian inference for stable distributions, Journal of the American Statistical Association, 90, 605-613. [8] Carr, P., H. Geman, D.B. Madan and M. Yor, 2002, The Fine structure of asset returns: An empirical investigation, Journal of Business, vol. 75, no.2, 305-332. [9] Carr, P. and D.B. Madan, 1999, Option valuation using the fast Fourier trans- form, Journal of Computational Finance, 3, 463-520. [10] Carr, P. and L. Wu, 2003, The nite moment log stable process and option pricing, The Journal of Finance, vol. LVIII, no.2, 753-778. [11] Carr, P. and L.Wu, 2004, Time-changed L evy processes and option pricing, Journal of Financial Economics, 71, 113-141. 90 [12] Chambers, J.M., C.L. Mallows and B.W. Stuck, 1976, A method for simulating stable random variables, Journal of the American Statistical Association, vol. 71, no. 354, 340-344. [13] Chernov, M. and E. Ghysels, 2000, A study towards a unied approach to the joint estimation of objective and risk neutral measures for the purpose of option valuation, Journal of Financial Economics, 56, 407-458. [14] Chib, S. and E. Greenberg, 1996, Markov Chain Monte Carlo simulation meth- ods in econometrics, Econometric Theory, 12, 409-431. [15] Chib, S., F. Nardari and N. Shephard, 2002, Markov chain Monte Carlo meth- ods for stochastic volatility models, Journal of Econometrics, 108, 281-316. [16] Cox, J.C., J.E. Ingersoll, Jr., and S.A. Ross, 1985, A theory of the term struc- ture of interest rates, Econometrica 53, 385-407. [17] Devroye, L., 1986, Nonuniform random variate generation (Springer Verlag, New York) [18] Durham, G., 2004, Monte Carlo methods for Estimating, Smoothing and Filter- ing One- and Two-Factor Stochastic Volatility Models, University of Colorado, working paper. [19] Engle, R.F., 1982, Autoregressive conditional heteroscedasticity with estimates of the variance of United Kingdom in ation, Econometrica 50, 987-1007. [20] Eraker, B., M. Johannes and N. Polson, 2003, The impact of jumps in volatility and returns, The Journal of Finance, vol. LVIII, no.3, 1269-1300. [21] Gallant, A.R. and G. Tauchen, 1996, Which moments to match?, Econometric Theory, 12, 657-681. [22] Geman, S. and D. Geman, 1984, Stochastic relaxation, Gibbs distributions and the Bayesian restoration of images, IEEE Transactions on Pattern Analysis and Machine Intelligence, 12, 609-628. 91 [23] Geweke, J., 1995, Monte carlo simulation and numerical integration, Sta re- port no. 192 (Federal Reserve Bank of Minneapolis, MN). [24] Gilks, W.R. and P. Wild, 1992, Adaptive rejection sampling for Gibbs sampling, Applied Statistics, vol. 41, no. 2, 337-348. [25] Hamilton, J.D., 1994, Time series analysis (Princeton University Press, Prince- ton, NJ). [26] Heston, S., 1993, A closed - form solution for options with stochastic volatility with applications to bonds and currency options, Review of Financial Studies, 6, 327-344. [27] Huang, J.Z. and L. Wu, 2004, Specication analysis of option pricing models based on time-changed L evy processes, Journal of Finance, vol. LIX, No.3, 1405-1439. [28] Jacquier, E., N.G. Polson and P.E. Rossi, 1994, Bayesian analysis of stochastic volatility models, Journal of Business & Economic Statistics, vol. 12, no. 4, 371-389. [29] Jacquier, E., N.G. Polson and P.E. Rossi, 2004, Bayesian analysis of stochastic volatility models with fat-tails and correlated errors, Journal of Econometrics, 122, 185-212. [30] Janicki, A. and A. Weron, 1994, Simulation and chaotic behavior of -stable stochastic processes (Marcel Dekker, New York) [31] Johannes, M. and N. Polson, 2003, MCMC methods for Financial Econometrics, 2003, Working Paper [32] Jones, C., 1998, Bayesian estimation of continuous - time Finance models, Working Paper [33] Jones, C.S., 2003, The dynamics of stochastic volatility: evidence from under- lying and options market, Journal of Econometrics 116, 181-224. 92 [34] Yu, J., 2005, On leverage in a stochastic volatility model, Journal of Econo- metrics, 127, 165-178. [35] Kloeden, P.E. and E. Platen, 1992, Numerical simulation of stochastic dier- ential equations (Springer-Verlag, New York, NY). [36] Li, H., M.T. Wells and C.L. Yu, 2006, A Bayesian Analysis of Return Dynamics with Stochastic Volatility and Levy Jumps, Working Paper [37] Lombardi, M., 2004, Bayesian inference for alpha-stable distributions: A ran- dom walk MCMC approach, University of Florence, Working Paper [38] Madan, D.B., P. Carr and E.C. Chang, 1998, The variance gamma process and option pricing, European Finance Review, 2, 79-105. [39] Merton,R.C.,1973, Theory of rational option pricing, Bell Journal of Economics and Management Science 4, 141-183. [40] Metropolis, N., A.W. Rosenbluth, M.N. Rosenbluth, A.H. Teller and E. Teller, 1953, Equations of state calculations by fast computing machines, Journal of Chemical Physics, 21, 1087-1092. [41] Nolan, J.P., 2005, Stable Distributions, models for heavy tailed data, American University Manuscript [42] Pitt, M. and N. Shephard, 1999, Filtering via Simulation: Auxiliary Particle Filter, Journal of the American Statistical Association, 94, 590-599. [43] Ripley, B.D., 1987, Stochastic simulation (John Wiley, New York) [44] Roberts, G. and J. Rosenthal, 1998, Markov Chain Monte Carlo: Some practical implications of theoretical results, Canadian Journal of Statistics, 26, 4-31. [45] Rubinstein, M., 1994, Implied binomial trees, Journal of Finance, 49, 771-818. [46] Samorodnitsky, G. and M.S. Taqqu, 1994, Stable non-Gaussian random pro- cesses: stochastic models with innite variance (Chapman & Hall, New York) 93 [47] Vasicek, O., 1977, An equilibrium characterization of the term structure, Jour- nal of Financial Economics 5, 177-188. [48] Zolotarev, V.M., 1986, One - dimensional stable distributions, Amer. Math. Soc. Transl. of Math. Monographs, vol. 65, American Mathematical Society, Providence, RI (translation of the original 1983 Russian) 94 Appendix A Derivation of complete conditional distributions Models with stochastic volatility via diusion only: SVPJ, SVPJLE, SVSJ, SVSJLE, SVVG, SVVGLE - non-jump specic parameters complete condi- tional posteriors. In the following, let denote = nfg for sets and , residuals from log- volatility equation " (3) t = " (1) t + p 1 2 " (2) t , residuals from returns equation " (1) t and J t is the jump size dened as: J t = 8 > > > > > > < > > > > > > : q t t for models with Poisson jumps J VG t for models with VG jumps S t (;; 0; j 1= ) for models with stable jumps In case of models without leverage eect we impose restriction on = 0. Let NJ = (; h ; h ; h ;) 95 be the vector of parameters (non-jump specic parameters), X NJ = (fh t g T t=1 ) be the vector of latent log-volatility states, Y = (Y 1 ;Y 2 ;:::;Y T ) =fY t g T t=1 be the vector of observed log-asset prices or index levels and J = (J 2 ;:::;J T ) =fJ t g T t=2 be the jump sizes of respective models. 1. updating : We choose the following prior distribution on : N(a;A): In our application we set a = 0 and A = 10, which is relatively at prior for the asset returns mean. The conditional posterior distribution is conjugate to prior and also normal given by: p(j h ; h ; h ;;X NJ ;Y;J)N(a ;A ); where A = (A 1 + 1 2 P T1 t=1 1 exp(ht) ) 1 and a =A (A 1 a + 1 1 2 P T1 t=1 h (Y t+1 YtJ t+1 ) expht " (3) t+1 0:5 exp(0:5ht) i : 96 2. updating h : The prior on h is h truncated N(b;B) such that h 2 (0; 2 ) with b = 1 and B = 6 in our application, which is also a relatively at prior that imposes stationarity on the log-volatility process h t . Hence, the conditional posterior is also truncated at 0 and 2 and is conjugate to prior: p( h j; h ; h ;;X NJ ;Y;J)N(b ;B ); and h 2 (0; 2 ); where B = (B 1 + 1 1 2 P T1 t=1 ( h ht) 2 2 h ) 1 and b =B (B 1 b + 1 2 h (1 2 ) P T1 t=1 h (h t+1 h t )( h h t ) h 1=2 " (1) t+1 ( h h t ) i . 3. updating h : The prior on h is h N(c;C); in our application c = 0 and C = 10, which is also a relatively at prior. Hence, the conditional posterior is conjugate to prior: p( h j; h ; h ;;X NJ ;Y;J)N(c ;C ); where C = (C 1 + (T1) 2 h 2 h (1 2 ) ) 1 and c =C (C 1 c + h 2 h (1 2 ) P T1 t=1 h h t+1 h t + h h t h 1=2 " (1) t+1 i . 97 4. updating h and . We update h and in a block following Jacquier, Polson and Rossi (2004) ([29]). Let dene the following 1-1 correspondence: h = h ; (A.1) ! h = 2 h (1 2 ): (A.2) We choose the joint prior distribution on the transformed parameters h and ! h specied by! h IG(d;D) and h j! h N(0; 1 2 ! h ) as in Jacquier, Polson and Rossi (2004) ([29]). In our application we choose d = 3 and D = 1 20 for an uninformative prior. The conditional posterior of ! h is conjugate to prior 1 and given by: p(! h j; h ; h ;X NJ ;Y;J)IG(d ;D ); where d = d + T1 2 ; D =D + 1 2 (err 0 err) +b 0 par ( P T1 t=1 (" (1) t+1 ) 2 + 2) 1 ( P T1 t=1 (" (1) t+1 ) 2 )b par ; whereerr is the vector of regression residuals and b par is an OLS estimator of h in the following regression model: h t+1 h t p h ( h h t ) p = h " (1) t+1 + p ! h t+1 ; t+1 N(0; 1); t = 1; 2;:::;T 1: 1 Note, that we do not condition on h : 98 The conditional posterior of h is also conjugate to prior distribution and is given by: p( h j! h ;; h ; h ;X NJ ;Y;J)N(m ;M ); where M =! h (2 + P T1 t=1 (" (1) t+1 ) 2 ) 1 , m =M P T1 t=1 (" (1) t+1 [ h t+1 ht p h ( h h t ) p ]): After draw of ( h ;! h ) we found ( h ;) by an inverse of the correspondence in eq. (A.1) and (A.2) given by 2 h =! h + 2 h , = h h : 5. updating volatility state h t for t = 1: By application of the Bayes rule, we have: p(h 1 j NJ ;X NJ(h 1 ) ;Y;J)/p(Y 2 j NJ ;X NJ ;Y (Y 2 ) ;J)p(h 1 jh 2 ; NJ ); where the second probability is given by the symmetry formula for autoregressive models (AR): p(h 1 jh 2 ; NJ )N( h h + (1 h )h 2 ; 2 h ); and the rst probability is given by: p(Y 2 j NJ ;X NJ ;Y (Y 2 ) ;J) / 1 p exph 1 exp 8 > < > : 1 2 Y 2 Y 1 J 2 (exph 1 ) 1 2 " (3) 2 2 (1 2 ) exph 1 9 > = > ; : We use the MH step to sample from this distribution with proposal density given by p(h 1 jh 2 ; NJ ). 99 6. updating volatility state h t for t =T . By application of the Bayes rule, we have: p(h T j NJ ;X NJ(h T ) ;Y;J)/p(Y T j NJ ;X NJ ;Y (Y T ) ;J)p(h T jh T1 ; NJ ); and after simplifying: p(h T j NJ ;X NJ(h T ) ;Y;J)N( h (exph T1 ) 1 2 ; 2 h (1 2 )); where =h T1 + h ( h h T1 ); =Y T +Y T1 + +J T : 7. updating volatility states h t for 1<t<T: By application of the Bayes rule, we have: p(h t j NJ ;X NJ(ht) ;Y;J)/ p(Y t+1 j NJ ;X NJ ;Y (Y t+1 ) ;J)p(Y t j NJ ;X NJ ;Y fYt;Y t+1 g ;J) p(h t+1 jh t ; NJ )p(h t jh t1 ; NJ ); where p(Y t+1 j NJ ;X NJ ;Y (Y t+1 ) ;J)/ 1 p exph t exp ( 1 2 [Y t+1 Y t J t+1 (exph t ) 1 2 " (3) t+1 ] 2 (1 2 ) exph t ) ; p(Y t j NJ ;X NJ ;Y fYt;Y t+1 g ;J)/ exp ( 1 2 [Y t Y t1 J t (exph t1 ) 1 2 " (3) t ] 2 (1 2 ) exph t1 ) ; 100 and p(h t+1 jh t ; NJ )p(h t jh t1 ; NJ )N(; 2 h (1 h ) 2 + 1 ); = (1 h )(h t+1 +h t1 ) + 2 h 2 h (1 h ) 2 + 1 : In this case we also use MH algorithm, with proposal density by p(h t+1 jh t ; NJ )p(h t jh t1 ; NJ ): We found this algorithm to be very ecient with acceptance rates of above 80% for most data sets. Jump specic parameters and latent variables of models with Poisson jumps SVPJ ( = 0) and SVPJLE. Let adopt similar notation to paragraph above. Moreover, dene J = f j ; j ; j g for jump specic parameters, X J =ffq t g T t=2 ;f t g T t=2 g for jump specic latent variables and X =X NJ [X J and = NJ [ J . 1. updating latent jump timesfq t g T t=2 . Since prior onq t+1 isBernoulli( j ), the conditional posterior is also Bernoulli: p(q t+1 = 1j;X q t+1 ;Y )/p(Y t+1 jY t ;X q t+1 ;;q t+1 = 1)p(q t+1 = 1j j ) and p(q t+1 = 0j;X q t+1 ;Y )/p(Y t+1 jY t ;X q t+1 ;;q t+1 = 0)p(q t+1 = 0j j ) 101 where the sum of probabilities is normalized to one. In the above the rst probability is given by: p(Y t+1 jY t ;;X)/ exp( 1 2 [Y t+1 Y t q t+1 t+1 ( exph t ) 1 2 " (3) t+1 ] 2 (1 2 ) exph t ) and second by Bernoulli( j ) prior distribution. 2. updating latent jump sizesf t g T t=2 . By application of the Bayes rule, we have: p( t+1 j;X t+1 ;Y )/p(Y t+1 jY t ;X t+1 ;; t+1 )p( t+1 j j ; j ) Ifq t+1 takes value zero, then the rst component is not a function of t+1 jump size and we are left with prior distribution p( t+1 j j ; j ) which is normal with mean j and variance 2 j . If q t+1 takes value one, then both components are the function of t+1 and we have: p( t+1 j;X t+1 ;Y )N( 1 2 B A ; (1 2 ) exph t 2 j A ) where A = 2 j + (1 2 ) exph t and B =2 2 j (Y t+1 Y t ( exph t ) 1 2 " (3) t+1 ) 2 j (1 2 ) exph t : 102 3. updating jump intensity parameter j . We specify beta prior distribution on j Beta(0:5; 0:5). Then by the Bayes rule the conditional posterior is also beta distributed: p( j j j ;X;Y )/p(fq t g T t=2 j j )p( j ) hence p( j j j ;X;Y )Beta( T X t=2 q t + 0:5; T X t=2 (1q t ) + 0:5): 4. updating mean jump size parameter j : By application of the Bayes rule, we have: p( j j j ;X;Y )/p(f t+1 gj j ; 2 j )p( j ): Using standard OLS theory we take normal prior with mean a = 0 and variance A = 10, which is relatively at prior, and have the following conditional posterior: p( j j j ;X;Y )N(a ;A ); where A = (A 1 + (T 1) 2 j ) 1 and a =A (A 1 a + 2 j T P t=2 t ). 5. updating variance jump size parameter 2 j : By application of the Bayes rule, we have: p( 2 j j 2 j ;X;Y )/p(f t+1 gj j ; 2 j )p( 2 j ): 103 Using standard OLS theory we take inverse gamma prior on 2 j IG(b;B) with b = 3 and B = 1 20 , which is relatively at prior, and have the following conditional posterior: p( 2 j j 2 j ;X;Y )IG(b ;B ); where b = (T 1)=2 +b and B =B + 0:5 T P t=2 ( t j ) 2 . Jump specic parameters and latent variables of models with Variance Gamma jumps SVVG ( = 0) and SVVGLE. In this paragraph we redene the jump specic vector of parame- ters J =( j ; j ; j ) and jump specic vector of latent variables X J = (fG t g T t=2 ;fJ VG t g T t=2 ). All other notation is left as in the previous part of the ap- pendix. 1. updating j : We impose relatively at normal prior on j with mean m = 0 and variance M = 100. By standard OLS theory complete conditional posterior is also normal and given by: p( j j j ;X;Y )N(m ;M ); where M = (M 1 + 2 j T P t=2 G t ) 1 and m =M (M 1 m + 2 j T P t=2 J VG t ). 104 2. updating 2 j : We impose relatively at prior on 2 j IG(c;C) with c = 3 and C = 1 20 . By standard OLS theory the complete conditional posterior is also inverse gamma: p( 2 j j 2 j ;X;Y )IG(c ;C ) with c = (T 1)=2 +c and C =C + 0:5 T P t=2 ( J VG t p Gt j p G t ) 2 . 3. updating j : We impose completely uninformative prior on j U(0; 1000) and truncate the prior for j > 1000 values. By application of the Bayes rule for 0 < j < 1000 we have: p( j j j ;X;Y )/p( j jfG t g T t=2 )/p(fG t g T t=2 j j )p( j )/p(fG t g T t=2 j j )/ T Y t=2 2 6 4 1 ( j ) j ( j ) G j t 3 7 5 exp[ 1 j T X t=2 G t ] We draw from this distribution by using MH random walk algorithm. 4. updatingfG t g T t=2 : By application of the Bayes rule, we have: p(G t j;Y;X Gt )/p(G t j;J t )/p(J t jG t ;)p(G t j); 105 where the rst component is normal and the second is gamma kernel. By substitu- tion of respective formulas, we have the following conditional posterior for G t : G ( j 3 2 ) t exp G t ( 1 j 2 j 2 2 j ) exp (J VG t ) 2 2 2 j G t : This density is unimodal but not log-concave for all range of parameters. We use MH random walk algorithm to sample from this distribution. 5. updatingfJ VG t g T t=2 : By application of the Bayes rule, we have: p(J VG t j;X J VG t ;Y )/p(Y t jY t1 ;J VG t ;fh t g;)p(J VG t j;G t ); where the both kernels are normal with respect to J t and hence the conditional posterior is also normal. The mean of the posterior is given by: 2 j G t C t + exph t1 (1 2 ) j G t 2 j G t + exph t1 (1 2 ) and variance by: 2 j G t exph t1 (1 2 ) 2 j G t + exph t1 (1 2 ) ; where C t =Y t Y t1 p exph t1 " (3) t . Updating parameters and latent variables in VGSV model specication with stochastic volatility from jumps only. 106 First note, that all non-jump specic parameters and latent variables can be updated as in model SJSV and this update procedure was already discussed. More- over update of j , 2 j is performed in the exactly same way as for models SVVG and SVVGLE above. In case of updatingfG t g T t=2 and j we have to replace value by exp( t1 ) in the respective formulas above, since we condition also on the value of log-volatility from jumpsf t g T t=2 . Hence, we only have to describe the update pro- cedure for log-volatiltiy from jumpsf t g T t=2 to completely characterize the MCMC algorithm. 1. updating log-volatiltiy from jumpsf t g T t=2 : First, for 0<t<T we have: p( t j;X t ;Y )/p( t j;fG t g; t+1 ; t1 )/ p(G t+1 j t ;)p( t j; t1 ; t+1 ); where the rst kernel is of gamma form and the latter of normal. Substituting all formulas, we have: G t+1 j exp t j 1 ( expt j ) exp ( t g) 2 2G ; whereg andG are respectively mean and variance parameters of the normal kernel part. Since in case of 1<t<T we can condition on both t+1 and t1 we have: g = (1 )( t+1 t1 ) + ( ) 2 (1 ) 2 + 1 107 G = 2 (1 ) 2 + 1 : However, for t = 1 we can condition only on 2 and hence we have: g = ( ) + (1 ) 2 G = 2 : Finally, fort =T we have much simpler conditional posterior, since we do not have gamma kernel in the formula: p( T j;X T ;Y )/p( T j; T1 ) and hence: p( T j;X T ;Y )N(g;G); where g = ( ) + (1 ) T1 G = 2 : For t<T we use MH random walk algorithm and for t =T we directly draw from normal distribution with parameters specied above.

Abstract (if available)

Abstract In this work we analyze asset returns models with diffusion part and jumps in returns with stochastic volatility either from diffusion or pure jump part. We consider different specifications for the pure jump part including compound Poisson, Variance Gamma and Levy alpha-stable jumps. Monte Carlo Markov chain algorithm is constructed to estimate models with latent Variance Gamma and Levy alpha-stable jumps. Our construction corrects for separability problems in the state space of the MCMC for Levy alpha-stable jumps. We apply our model specifications for analysis of S&P500 daily returns. We find, that models with infinite activity jumps and stochastic volatility from diffusion perform well in capturing S&P500 returns characteristics. Models with stochastic volatility from jumps cannot represent excess kurtosis and tails of returns distributions. One-day and one-week ahead prediction and VaR performance characterizing conditional returns distribution rejects Variance Gamma jumps in favor of Levy alpha-stable jumps in returns. JEL classification: C1

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Stochastic models: simulation and heavy traffic analysis

PDF

A general equilibrium model for exchange rates and asset prices in an economy subject to jump-diffusion uncertainty

PDF

Model, identification & analysis of complex stochastic systems: applications in stochastic partial differential equations and multiscale mechanics

PDF

Forward-backward stochastic differential equations with discontinuous coefficient and regime switching term structure model

PDF

Numerical analysis on the influence of information delay on optimal investment strategy in the family of 4/2 stochastic volatility models

PDF

Approximating stationary long memory processes by an AR model with application to foreign exchange rate

PDF

Three essays on linear and non-linear econometric dependencies

PDF

Essays on high-dimensional econometric models

PDF

Efficient inverse analysis with dynamic and stochastic reductions for large-scale models of multi-component systems

Asset Metadata

Creator Szerszen, Pawel (author)

Core Title Bayesian analysis of stochastic volatility models with Levy jumps

School College of Letters, Arts and Sciences

Degree Doctor of Philosophy

Degree Program Economics

Publication Date 04/21/2010

Defense Date 03/25/2008

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag asset returns,Bayesian estimation,Levy jumps,OAI-PMH Harvest,stochastic volatility

Language English

Advisor Magill, Michael J.P. (committee chair), Jones, Christopher S. (committee member), Moon, Hyungsik Roger (committee member)

Creator Email szerszen@usc.edu

Permanent Link (DOI) https://doi.org/10.25549/usctheses-m1170

Unique identifier UC1195935

Identifier etd-Szerszen-20080421 (filename),usctheses-m40 (legacy collection record id),usctheses-c127-61608 (legacy record id),usctheses-m1170 (legacy record id)

Legacy Identifier etd-Szerszen-20080421.pdf

Dmrecord 61608

Document Type Dissertation

Rights Szerszen, Pawel

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Repository Name Libraries, University of Southern California

Repository Location Los Angeles, California

Repository Email cisadmin@lib.usc.edu