Access the full text.
Sign up today, get DeepDyve free for 14 days.
Loading next page...
/lp/arxiv-cornell-university/modeling-microstructure-price-dynamics-with-symmetric-hawkes-and-5aEzu0QNJc
- ISSN
- 0165-1889
- eISSN
- ARCH-3346
- DOI
- 10.1016/j.jedc.2017.04.004
- Publisher site
- See Article on Publisher Site
Abstract
This study examine the theoretical and empirical perspectives of the symmetric Hawkes model of the price tick structure. Combined with the maximum likelihood estimation, the model provides a proper method of volatility estimation specialized in ultra-high-frequency analysis. Empirical studies based on the model using the ultra-high-frequency data of stocks in the S&P 500 are performed. The performance of the volatility measure, intraday estimation, and the dynamics of the parameters are discussed. A new approach of di usion analogy to the symmetric Hawkes model is proposed with the distributional properties very close to the Hawkes model. As a di usion process, the model provides more analytical simplicity when computing the variance formula, incorporating skewness and examining the probabilistic property. An estimation of the di usion model is performed using the simulated maximum likelihood method and shows similar patterns to the Hawkes model. 1 Introduction The extensive observations and analysis of ultra-high-frequency nancial data has become increasingly available due to the development of computing schemes, massive storage devices, and electronic trade systems in the nancial markets. The ultra-high-frequency data includes the price dynamics and various types of trade orders recorded in seconds or with a shorter time resolution. Therefore, there has been growing attention in the necessity for proper analysis and modeling of ultra-high-frequency nancial data among practitioners and theorists. One of the important subjects of modeling ultra-high-frequency data is the price dynamics in micro level with tick structures. To describe the micro structure of the price dynamics and order ows, the Hawkes process (Hawkes, 1971a,b) has been used to consider the non-time-homogeneous features of the duration between price changes or orders such as clustering and mutual e ect. The Hawkes process belongs to the class of point processes and is de ned by constructing the conditional intensity processes as a function of previous events. Hewlett (2006) examined the model of the arrival times of trades and the price impacts based on a symmetric bivariate Hawkes process. Large (2007) examined the market resilience after large trades using the limit order book data and mutually excited multivariate Hawkes processes. Bowsher (2007) introduced a generalized Hawkes model to analyze the relationship between the trading times and mid price changes. With mutually excited Hawkes processes that have a strong microscopic mean reversion property, Bacry et al. (2013) constructed a model that accounts for the market microstructure noise and the Epps e ect. On the other hand, Da Fonseca and Zaatour (2014b) focused on the clustering behaviors of trades using self-excited Hawkes processes with an application to the generalized method of moments estima- tions. Da Fonseca and Zaatour (2014a) provided the moment conditions and autocorrelation functions of self and mutually excited Hawkes processes to exhibit both clustering and mean reversion. Bacry and Muzy (2014) proposed a multivariate Hawkes process to model the price dynamics and the market impact of market orders to account for the various stylized facts of the market microstructure. For more previous nancial studies on market microstructure or price dynamics based on point processes or intensity modeling, the reader should refer to Bauwens and Hautsch (2009), Embrechts et al. (2011), Bacry et al. (2012), Zheng et al. (2014) and Choe and Lee (2014a). The Hawkes process has also been Department of Statistics, Yeungnam University, Gyeongsan, Gyeongbuk 38541, Korea Corresponding author, School of Management Engineering, UNIST(Ulsan National Institute of Science and Technol- ogy), Ulsan 44919, Korea arXiv:1908.05089v1 [q-fin.ST] 14 Aug 2019 applied to modeling the credit and contagion risk, see Errais et al. (2010), A t-Sahalia et al. (2010) and Dassios and Zhao (2012). This paper focuses on the tick price dynamics and volatility estimation. The realized volatility estimator (Barndor -Nielsen and Shephard, 2002a,b; Andersen et al., 2003) in the ultra-high-frequency dynamics can be biased; when one uses every sample of ultra-high-frequency nancial data to calculate the nite sum approximation of the integrated volatility due to the microstructure noise and clustering property, see Hansen and Lunde (2006). The adjustment methods of the bias in a nonparametric fashion (Zhang et al., 2005; A t-Sahalia et al., 2005, 2011) have been introduced. In these approaches, one supposes that the observed price process consists of the latent ecient price and noise term around the ecient price process. In contrast, in the Hawkes models or di usion approach introduced in this paper and related literatures, one models the observed price movements directly, which may include the noise, and compute the closed form formula for the variance of the return and analyze the properties of the variance. Empirical studies to compare the volatilities calculated by Hawkes modeling and realized quadratic variation using the stock prices of the S&P 500 were performed. Because the Hawkes model approach incorporates all of the arrival times of the price change within a millisecond time resolution, such richness of data provides the eciency of the volatility estimation. This paper reports the relative eciency of the Hawkes volatility compared to the realized volatility in simulation studies. Therefore, owing to the rich information in ultra-high-frequency data combined with ecient likelihood estimation methods, one can estimate the parameters and volatilities within a relatively short time period of observation. This is one of important features of the Hawkes model, and with this property, this paper presents the empirical results of the intraday volatility dynamics based on the Hawkes model. By observing the intraday volatility variation in every moment, one can respond to sudden market movements more e ectively. In addition, a di usion counterpart of the Hawkes model for the micro price dynamics is introduced. The di usion model consists of the square root processes for both volatility and drift. The proposed di usion model has similar properties to the symmetric Hawkes model of price process such as the strong correlation of the mean process over the time lag on a small time scale and hence it incorporates the market microstructure noise. This paper reports that the di usion models generate the distribution very close to the corresponding Hawkes models using the Kolmogorov forward equation. As a di usion model, it is simpler to compute variance formula, able to introduce the leverage parameter which explain the skewness and provides the insight about the distributional property of return. In addition, using simulated likelihood estimation method, the model parameters and volatility of the equity returns are examined. The remainder of the paper is organized as follows. Section 2 introduces the Hawkes model for the micro price dynamics with the basic setup similar to Hewlett (2006). Section 3 proposes and discusses the di usion analogy of the symmetric Hawkes model. Section 4 shows the empirical results with the symmetric Hawkes model and the corresponding di usion model. The daily and intraday variation of the Hawkes parameters and volatility with several stock data of the S&P 500 are shown. Section 5 concludes the paper. The proofs and further explanations are gathered in the Appendix. 2 Hawkes process for tick dynamics 2.1 Point process This section starts with the introduction of the Hawkes process, which belongs to the class of point processes, (see, Daley and Vere-Jones (2003)). A point process, N , is formally de ned on a state space, X , as a mapping from a probability space ( ;P) to N , where N denotes the space of all counting measures on the - eld of X 's Borel sets, B . The space X is a complete separable metric space and to study the tick-dynamics of a stock price movements, this paper focuses on the case that X = R, the time domain. As a counting measure, N (A; !) has a non-negative integer value for any measurable set A 2 B and is nite for any bounded measurable A. Using the Dirac measure, , de ned for every X x x 2 X , the counting measure is represented by N = k i x 2 where fx g is a countable set with at most nitely many x in any bounded Borel set and k is a positive i i i integer. This paper only considers the simple counting measure, i.e., k = 1 for all i. A point process N can be regarded as a stochastic process by letting N (t; !) = N (( 1; t]; !). Consider a ltered probability space ( ;fF g;P), 1 < t T , where the - eld F is generated t t by N (t). The Hawkes process is an orderly stationary point process N constructed by modeling the conditional intensity, . The conditional intensity function is represented as an adapted process to fF g such that (t)dt = E[N (t + dt) N (t)jF ]. For an M -dimensional Hawkes process (N ; : : : ; N ), each t 1 M intensity, (t) of N is assumed to be i i (t) = + (t u)dN (u) i i i;j j j=1 where (t u) is normally a deterministic function and called kernel. The integration of the r.h.s. is i;j the stochastic integration de ned pathwise. To apply the stochastic integration theory in the later, the Hawkes and intensity processes are considered to be right continuous processes with left limits. 2.2 Self and mutually excited Hawkes This subsection brie y reviews the self and mutually excited Hawkes model. Consider a two dimensional Hawkes process (N ; N ) with exponential decay kernels in the conditional intensities with constants , 1 2 i and , for 0 < t: ij ij Z Z t t (t u) (t u) 11 12 (t) = + e dN (u) + e dN (u) 1 1 11 1 12 2 1 1 Z Z t t t t (t u) (t u) 11 12 11 12 = + (0)e + (0)e + e dN (u) + e dN (u); (1) 1 11 12 11 1 12 2 0 0 and Z Z t t (t u) (t u) 21 22 (t) = + e dN (u) + e dN (u) 2 2 21 1 22 2 1 1 Z Z t t t t (t u) (t u) 21 22 21 22 = + (0)e + (0)e + e dN (u) + e dN (u) (2) 2 21 22 21 1 22 2 0 0 where (t u) ij (t) = e dN (u): ij ij j In this paper, this model is called the fully characterized self and mutually excited Hawkes process compared to the symmetric Haweks process introduced later. Note that and are self-excited 11 22 components, and are mutually excited components, and every parameter such as and , 12 21 ij ij can have a di erent value. This model was proposed by Bacry et al. (2013) and was studied for a simpli ed version focused on the self-excited term. The self and mutually excited Hawkes model and its moment properties are studied in Da Fonseca and Zaatour (2014a). The components of the intensity processes, , can be rewritten by ij (t u) ij (t) = q e dN (u) (3) ij ij ij j ij (t u) ij where q := and the integrand, e , is a normalized decaying function in the sense that ij ij ij ij e d = 1: ij The coecients, q , form a branching matrix, Q = fq g and if the spectral radius, the maximum ij ij i;j=1;2 of the absolute eigenvalues of Q, is less than 1, then the Hawkes process is well de ned (Hawkes and Oakes, 1974; Br emaud, 1981). 3 The stock price process can be assumed to be represented by the di erence between two Hawkes processes, S = S + fN (t) N (t) (N (0) N (0))g (4) t 0 1 2 1 2 where denotes the unit size of the price movement in the tick structure of price dynamics. (In the previous subsection, was used to denote the Dirac measure. On the other hand, without the subscript, is a constant that represents the tick size.) The process N represents the up movements of the price process and N represents the down movements. However, the fully characterized Hawkes model is too complicated not only in the number of pa- rameters but also in the fact that the model becomes four dimensional problems when dealing with the moment conditions as explained in A. (Nonetheless, we will provide some empirical results with the fully characterized model in Section 4.) In the next subsection, we consider a simpler version. 2.3 Symmetric Hawkes process This subsection explains the symmetric Hawkes model for the price dynamics. The empirical study shows that the symmetric version also well represents the basic properties of the tick dynamics. To simplify the model from the fully characterized version, the parameter condition is imposed as := = ; := = c 12 21 s 11 22 := = = = ; := = : 11 12 21 22 1 2 Then Z Z t t (t u) (t u) (t) = + e dN (u) + e dN (u) (5) 1 s 1 c 2 1 1 Z Z t t t (t u) (t u) = + ( (0) )e + e dN (u) + e dN (u) 1 s 1 c 2 0 0 Z Z t t (t u) (t u) (t) = + e dN (u) + e dN (u) (6) 2 c 1 s 2 1 1 Z Z t t t (t u) (t u) = + ( (0) )e + e dN (u) + e dN (u): 2 c 1 s 2 0 0 This can also be written as d (t) = ( (t))dt + dN (t) + dN (t) 1 1 s 1 c 2 = f + ( ) (t) + (t)g dt + (dN (t) (t)dt) + (dN (t) (t)dt) s 1 c 2 s 1 1 c 1 1 d (t) = ( (t))dt + dN (t) + dN (t) 2 2 c 1 s 2 = f + (t) + ( ) (t)g dt + (dN (t) (t)dt) + (dN (t) (t)dt): c 1 s 2 c 1 1 s 1 1 Note that (t) = (t); (t) = (t): c 11 s 21 s 12 c 22 By setting = and = , the processes (N ; N ; ; ) are Markov and the di erential 11 12 21 22 1 2 1 2 equation system of the expected intensities becomes two dimensional. By the di erential forms of , ` (tjs) ` (tjs) s c 1 = + (7) ` (tjs) ` (tjs) c s 2 where ` (tjs) = E [ (t)] and the derivatives are with respect to t. Let i s i s c M = : c s The eigenvalues of M are ( ; ) = ( + ; + + ); 1 2 c s c s 4 and the corresponding eigenvectors are ( 1; 1) and (1; 1), respectively. If the eigenvalues are all negative, then the solution to the system converges to the particular solution as time approaches in nity. This is equivalent to the condition that the spectral radius of the branching matrix is less than one where, in the sense of parametrization in Eq. (3), the branching matrix is q q s c Q = q q c s with q := = and q := = . s s c c The solution of system (7) is (s) + (s) (s) + (s) E [ (t)] 1 1 1 1 2 1 2 s 1 (t s) (t s) (t s) 1 2 2 = e + e 1 e : E [ (t)] 1 1 1 s 2 2 2 The long-run expectations of the intensities as t ! 1, i.e., the particular solution of the system (7) is 1 1 = : 1 1 ( + ) s c 2 In the latter, for computational ease, it is usually assumed that the intensity processes are in the sta- tionary state at time 0, i.e., (0) = (0) = = : (8) 1 2 s c 2 The formula for the variance of the return generated by the symmetric Hawkes model is quite simple, as represented in Proposition 3. The simplicity largely depends on the symmetry of the parameter setting and the assumption of the stationary state condition at time 0. Indeed, the stationary condition does not signi cantly a ect the result on the variance formula in the high-frequency price dynamics modeling as the expectations of the intensities quickly converge. The formula was derived independently but Da Fonseca and Zaatour (2014a) reported a similar result (however, di erent from the exponential term below). Proposition 1. Assume that the price process, S, follows the di erence of two symmetric Hawkes processes de ned by Eqs. (4),(5), and (6). Under the stationarity condition of the intensity processes at time 0 as in Eq. (8), the variance of the return is represented by 2 t 2 t 1 1 S S 2 (0) e 1 e 1 t 0 1 2 2 Var = t 2( ) + ( ) s c s c 2 2 S S 2 0 1 1 0 1 Proof. See D. Remark 2. If t is suciently large, then the variance is approximated by 2 2 2 S S 2 (0)t 2 t t 0 1 Var = 2 2 2 S S S (1 q + q ) (1 q q ) 0 s c s c 0 1 0 2 t = : 2 s c s c S 1 + 1 In this approximation, the following parameterization of the symmetric Hawkes process is useful for the volatility estimation. The parameter is represented by a formula consisting of the annualized daily volatility and the other parameters in the Hawkes model and are given by 2 2 ann 1 3 2 T 2 2 2 (1 q + q ) (1 q q ) s c s c ann T 2 where denotes the annualized volatility, = =S and T is one year. ann r 0 5 Table 1: Simulation study with 500 samples H. vol TSRV s c True 0.0100 0.4000 0.5000 1.5000 0.1171 0.1171 mean 0.0100 0.4021 0.5027 1.5024 0.1177 0.1165 std. (0.0005) (0.0394) (0.0428) (0.0841) (0.0057) (0.0114) True 0.0500 0.6500 0.2000 1.7000 0.3396 0.3396 mean 0.0500 0.6514 0.2012 1.7027 0.3400 0.3370 std. (0.0014) (0.0282) (0.0144) (0.0643) (0.0103) (0.0283) Remark 3. Proposition 2 in Da Fonseca and Zaatour (2014a) showed the formula for the mean signature plot: (1 e 2 2 + (1 ) : Based on the de nition of the mean signature plot, by setting = t and multiplying the mean signature S S 2 t 0 plot by t=S(0) , the meaning of the formula is the same as the Var in Proposition 1 of our paper. If we rewrite Da Fonseca and Zaatour (2014a)'s formula with the notations used in our paper, then we have 2 t 2 (0) e 1 2 2 t + ( ) 2 ( ) : s c s c 2 2 0 1 The result of this formula is di erent from the formula for Proposition in the exponential term. However, as discussed in Remark 2, the exponential term is negligible if t is large. 2.4 Simulation study In this subsection, simulation studies are performed with the symmetric Hawkes processes. With prede- termined parameter settings, 500 sample paths of the price processes de ned by the di erence between the two symmetric Hawkes processes with 5.5 hours' time horizon are generated. For each path, the max- imum likelihood estimation is performed using the realized arrival times of the simulated path. Table 1 lists the results. The detailed information about the simulation method, see B and for the likelihood estimation, see C. The table consists of two panels with di erent parameter settings. The row `mean' is for the sample mean of the likelihood estimates of 500 samples. The row `std.' is for the sample standard deviations of the estimates. The column `H. vol' is for the mean of the volatility estimates calculated by the likelihood estimates of ; ; ; using Proposition 3. This is compared with s c the theoretic volatility computed by Proposition 3 in the row of `True'. The column `TSRV' reports the two scale realized volatility (TSRV) proposed by Zhang et al. (2005), which is known to be an unbiased estimator in the presence of independent market microstructure noise. For the TSRV computation, the small time scale is 1 second and the large time scale is 5 minutes. The Hawkes volatility and TSRV both are quite close to the true value of the volatility. The standard deviations of the Hawkes volatility are smaller than the standard deviations of the TSRV, implying the eciency of the maximum likelihood estimation. More precisely, in the maximum likelihood estimation of the Hawkes model, all the information about the time arrivals of events are used without missing single events over the observed period. On the other hand, in the computation of the realized volatility under the equidistant setting, it is needed to choose speci c points that belong to the sub-grids of the interval. The likelihood function of the symmetric Hawkes model may not be concave but is concave when is xed. For any given observed jump times t , the log likelihood function of the up jump over interval 6 [0; T ] is Z Z T T log L (T ) = log (u)dN (u) (u)du 1 1 1 1 0 0 Z Z t T = log (t ) (u)du (u)du 1 i 1 1 t t i 1 t <T X (T t ) i N e 1 e 1 = log (t ) (t ) (T ) 1 i 1 i 1 t <T e 1 where t is the last jump time up to T and = t t . Using Eq. (5), the term log (t ) (t ) N i i+1 i 1 i 1 i is represented by e 1 log (t ) (t ) 1 i 1 i Z Z t t i i t t (t u) (t u) i i i i = log (0)e + (1 e ) + e dN (u) + e dN (u) 1 s 1 c 2 0 0 Z Z t t i i e 1 t t (t u) (t u) i i i i (0)e + (1 e ) + e dN (u) + e dN (u) : 1 s 1 c 2 0 0 When is xed, then the term is represented by i i 1 e e 1 log (t ) (t ) = log(c + c + c + c ) (c + c + c + c ) 1 i 1 i i;0 i;1 i;2 s i;3 c i;0 i;1 i;2 s i;3 c for some constants c ; c ; c ; and c . By simple calculation, we have the negative semide nite Hessian i;0 i;1 i;2 i;3 matrix of the term with respect to ; ; is s c 2 3 c c c c c i;1 i;2 i;1 i;3 i;1 4 5 H = c c c c c : i;1 i;2 i;2 i;3 1;i i;2 (t ) i 2 c c c c c i;1 i;3 i;2 i;3 i;3 Similarly, we de ne H , and the Hessian matrix of the log L(T ) is H = (H +H ) which is also 2;i 1;i 2;i t <T negative semide nite and implies the log-likelihood function is conditionally concave when is xed. This means that if we compute the log-likelihood for every value of a reasonable set of , (a numerical procedure will perform this task well, because the log-likelihood function is concave for any xed ), and by comparing the computed values, we can nd the maximum log-likelihood. Therefore, we set a possible interval for , for example, 2 [1; 3], and with suciently small step size, for example, 0.0001, we can nd the estimates which make the log-likelihood close enough to the maximum log-likelihood. The examples with the above simulation set of Table 1 is shown in Figure 1. Although the above method can guarantee nding the maximum value, because it is time-consuming, we generally use the numerical procedure such as the BFGS algorithm based on the newton method to nd the maximum likelihood estimates. Although the proof of the BFGS algorithm's global convergence for the nonconvex function is not yet known, it is also known that global convergence works well in most cases (Li and Fukushima, 2001). For more information about the algorithm and its implementation, consult Broyden (1970) and Nash et al. (2014). In this simulation study, two methods show quite close results. For simulation set 1, the estimates computed by xing are = 0:0099; = 0:6590; = 0:0:4864; = 2:0646 and the estimates through s c the BFGS algorithm are = 0:0099; = 0:6590; = 0:4864; = 2:0346. For simulation set 2, the s c estimates computed by xing are = 0:0502; = 0:6273; = 0:2085; = 1:6861 and the estimates s c through the BFGS algorithm are = 0:0502; = 0:6272; = 0:2084; = 1:6860. In other simulation s c examples not recorded here, the BFGS algorithm always yields very similar results when compared with the method of xing . Since the method of xing is relatively time-consuming, by assuming that the BFGS algorithm provides very accurate estimates, we use the BFGS algorithm in future estimations. 7 4 ×10 -3000 -1.085 -1.09 -3040 -1.095 -3080 1 2 3 1 2 3 β β Figure 1: Maximum log-likelihood function when is xed for simulation set 1 (left) and 2 (right) 3 Di usion analogy 3.1 Di usion model This subsection proposes a new di usion approach for the tick structure. The di usion model is analogous to the symmetric Hawkes model and has a similar probabilistic property. When the price process is represented by the di erence of the two Hawkes process, the increment of the price process can be rewritten as S(t) = f(N (t) N (t))g 1 2 p p p p N (t) (t)t N (t) (t)t 1 1 2 2 = ( (t) (t))t + (t) p t (t) p t: 1 2 1 2 (t)t (t)t 1 2 Based on the empirical studies, a sucient number of price changes were observed during, e.g., one minute, and hence the normal approximation to the Poisson distribution N (t) (t)t i i N (0; 1) (t)t can be considered. Therefore, it is natural to consider a di usion analogy to the symmetric Hawkes model such as p p p p dN (t) (t)dt dN (t) (t)dt 1 1 2 2 p p (t) (t) (t)dB (t) + (t)dB (t) 1 2 1 1 2 2 (t) (t) 1 2 for some independent Brownian motions B and B . 1 2 By the independence, the in nitesimal variance is p p Var (t)dB (t) + (t)dB (t) = ( (t) + (t))dt t 1 1 1 2 1 2 which can be written p p p (t) + (t)dW = (t)dB (t) + (t)dB (t) 1 2 1 1 2 2 s 2 for some Brownian motion W . In the left hand side, ( (t) + (t)) as the instantaneous variance V 1 2 t of the price process. In addition, by treating ( (t) (t)) as the mean process n of the price process, 1 2 t a di usion analogy of the price process can be derived as follows: dS = n dt + V dW : (9) t t t Now the di usion analogies of n and V are constructed. This is because, by the de nition of (t), t t i df( (t) (t))g = ( )( (t) (t))dt 1 2 s c 1 2 ( ) p p dN (t) (t)dt dN (t) (t)dt 1 1 2 2 + ( ) (t) p (t) p ; s c 1 2 (t) (t) 1 2 log-likelihood log-likelihood let dn = (a a b)n dt + (a a ) V dW t s c t s c t =: n dt + V dW 1 t t where a ; a ; b are the di usion counterparts of ; ; , respectively. s c s c The micro structure of the price dynamics are slightly di erent from the macro dynamics as the non-zero drift term in the price process is observed. The drift term in the micro dynamics is also called the microstructure noise and related to the mutually excited feature in the Hawkes model. In addition, because 2 2 2 df ( (t) + (t))g = f2 + ( + ) ( (t) + (t))gdt 1 2 s c 1 2 ( ) p p dN (t) (t)dt dN (t) (t)dt 1 1 2 2 p p + ( + ) (t) + (t) ; s c 1 2 (t) (t) 1 2 under the similar argument of the drift, such as V = ( (t) + (t)); t 1 2 p p p V dW = (t)dB (t) (t)dB (t) t 1 1 2 2 p p dN (t) (t)dt dN (t) (t)dt 1 1 2 2 = (t) p + (t) p ; 1 2 (t) (t) 1 2 the familiar square root variance process as introduced in Heston (1993) is derived: 2bm dV = (b a a ) V dt + (a + a ) V dW t s c t s c t b a a s c =: ( V )dt + V dW : 2 t t where m is the di usion counterpart of . s v In this reasoning, the correlation that satis es d[W ; W ] = dt is represented by t t (t) (t) 1 2 = : (t) + (t) 1 2 In addition, if there is no jump for a suciently long interval and hence (t) ! and (t) ! , 1 1 2 1 then ! 0. Even though is represented by s or converges to zero, we consider that this constraint is better to be relaxed for exibility of the model. Note that the above derivation is not an exact mathematical justi cation, but rather to provide an intuition to construct a di usion model for the micro structure of price dynamics. For example, if necessary, the asymmetry in the price dynamics is simply introduced by a constant leverage parameter such that s v d[W ; W ] = dt as in the typical macro level price dynamics modeling. Overall, the price, mean and variance process are as follows: dS = n dt + V dW ; t t t dn = n dt + V dW ; t 1 t t t v s v dV = ( V )dt + V dW ; d[W ; W ] = dt t 2 t t t with the parameter relations: = b a + a 1 s c = b a a 2 s c 2bm b a a s c = (a + a ) s c = a a = + : s c 1 2 9 0.8 0.4 997 998 999 1000 1001 1002 1003 Figure 2: Numerically computed probability density function of the price driven by the di usion model and a histogram of the Hawkes model price by the simulation with 30 seconds (right) Note that with this analogy, corresponds to in the symmetric Hawkes model and corresponds 1 1 2 to . The di usion model is not an exact mathematical limit of the (symmetric) Hawkes model but has a very close distributional property with the Hawkes model. For recent studies about the limit theorem of the Hawkes process, consult Jaisson et al. (2015). 3.2 Basic property The di usion model has several advantages. First, by the forward Kolmogorov equation, the joint probability density function f (s; n; v; t) of the di usion model with s = S ; n = n ; v = V at time t t t t satis es the following partial di erential equation 2 2 2 2 2 @f @f @ @ v @ f v @ f @ = n + nf ( v)f + + + vf 1 2 2 2 2 @t @s @n @v 2 @s 2 @n 2 @v 2 2 2 @ f @ @ + v + vf + vf @s@n @n@v @s@v and the density function can be computed via a numerical procedure such as nite di erence method. More precisely, because the variable s comes up only in the derivative operators in the above equation, to reduce the dimension of the PDE, consider the Fourier transform of f with respect to s. That is i s f (n; v; t; ) = f (s; n; v; t)e ds @f @ f ^ ^ and the Fourier transforms of and are i f and f , respectively. Thus, by applying the Fourier @s @s transform to the PDE, ^ ^ 2 2 ^ ^ 2 2 ^ @f @f v @ f @f @ f = ( + n + i v) + + ( v) + + i v + v 1 2 2 2 @t @n 2 @n @v 2 @v 2 2 @ f + v + i ( n) v + + f: 1 2 @n@v 2 The transformed function f can be computed by numerical procedures and the probability density func- tion of the price is generated by applying inverse Fourier transform to the computed f . The distributions of the di usion model and the simulated histograms of the corresponding Hawkes process are compared in Figure 2. The parameter settings are = m = 0:09; = a = 0:6; = a = 0:3; = b = 2:5; = 0:2; S = 1000 s s c c 0 and the time horizon is 30 seconds. Second, the derivation of the variance formula in the di usion model is relatively simple compared to the symmetric Hawkes model due to the analytical simplicity of the di usion processes. To derive 10 Diffusion model Hawkes model 0.012 0.012 0.006 0.006 0 0 10:00 12:00 14:00 15:30 10:00 12:00 14:00 15:30 (a) Volatility by the di usion model (b) Volatility by the symmetric Hawkes model Figure 3: The comparison between the volatility computed by the di usion model and the symmetric Hawkes model the variance formula of the return, for simplicity, it is assumed that the variance process V is in the stationary state at time 0. This is a similar assumption that the intensity processes in the symmetric Hawkes model are in the stationary state at time 0. If V is in the stationarity state at time 0, then 2bm V = = b a a s c and since E[V ] = , 2bm t E[V ]ds = : b a a s c Similarly, if the mean process n is in the stationarity state at time 0, then n = 0. t 0 Proposition 4. Assume that the price process S follows Eq. (9) and the instantaneous variance V and mean processes n are in the stationary state at time 0. The variance of the return is ( ) 2 2 t t t 1 1 1 e + 4e 3 + 2 t S S 1 2 ( t 1 + e ) t 0 1 1 Var = + + t : 2 3 2 S S 2 0 1 1 Proof. See E. If t is suciently large, then the variance is approximated by 2 2 S S 1 t 2 t b t t 0 Var + + t = 2 2 2 2 S S S 0 1 0 1 0 1 which is analogous to Remark 2. When = 0, i.e., the drift of the price process is zero, the variance of the return is simply t=S . The volatility of the di usion model computed by Proposition 4 and the volatility of the symmetric Hawkes model computed by Proposition 3 were compared with the following parameter settings = a = 1:2; = a = 0:3; = b = 2:2; = m = 0:01; =S = 0:002 s s c c 0 in Figure 3. The volatilities were not annualized to show the increasing shape with time. The two volatilities are quite close to each other. Figure 4 shows the annualized volatility surface as a function of and with a xed = 4 10 . With a xed , with increasing = a a , (when the self-excited coecient a is larger than the 1 s c s mutually excited coecient a ) the volatility increases. This result is expected because the self-excited coecient is related to trade clustering. The increasing rate of the volatility with respect to depends on the level of . Because = b , 1 1 with a xed , a large implies a large b and a short persistence. In addition, a small implies a 1 1 small b and a long persistence. Therefore, when < 0, i.e., the mutually excited e ect is larger than the self excited e ect, a longer persistence (smaller ) of the mutually excited e ect implies a smaller volatility and a shorter persistence implies a larger volatility. On the other hand, when > 0, i.e., the self excited e ect is larger than the mutually excited e ect, a longer persistence of the self excited e ect Volatility Volatility 0.25 0.2 0.15 0.1 0.05 4 0.5 1 −0.5 Figure 4: Annualized volatility surface as a function of and implies a larger volatility and a shorter persistence implies a smaller volatility. This contrast is visualized in Figure 4 with the di erent sign of the slope of the volatility with respect to depending on whether > 0. Because there are many studies focused on the behavior of the realized variance in the presence of the market microstructure noise, this study also examined the realized variance under the di usion model. Consider a discretized time interval [0; T ] with time step . For convenience, let T= be an integer. The realized variance is de ned by the nite sum approximation to the quadratic variation of the return over a time interval. The signature plot over a xed interval is the realized variance over the interval de ned as a function of : T= C ( ) = (R R ) (n+1) n n=0 where R = (S S )=S is the return process. (Depending on the context, R could be the log-return t 0 0 process.) The above formula is indeed the de nition of the realized variance, which is the consistent estimator of the true variance of the return in the absence of microstructure noise by the semimartingale theory. On the other hand, empirical studies showed that the realized variance depends on the size of the partition due to the microstructure noise or clustering property (Hansen and Lunde, 2006; Da Fonseca and Zaatour, 2014b). For the di usion model, under the stationarity assumption of the time series of the squared return, (R R ) , and the mean signature plot is (n+1) C ( ) = E[C ( )] = E[(R R ) ] (n+1) ( ) 2 2 1 1 1 e + 4e 3 + 2 2 ( 1 + e ) = + + : 2 3 2 S 2 0 1 1 Figure 5 shows the mean signature plots with various in the left and in the right. For parameter settings, S = 1; = 2 10 and = 0:5 in the left and = 0:3 in the right. With negative values 0 1 of , implying a < a and a more pronounced self-excited e ect, the mean signature plot increases s c as approaches zero. On the other hand, with positive , implying a < a and a more pronounced s c self-excited e ect, the mean signature plot decreases as approaches zero. In both cases, when is too small, there is bias between the realized variance and the true variance, which is in contrast to the traditional understanding in statistics that a more exact result is obtained with a large sample size. In addition, with a suciently large , the expected realized variances converge. Third, as mentioned before, the asymmetry in the price distribution can be introduced easily with the leverage parameter . This method is a natural extension of the method used to introduce asymmetry in a macro level price dynamics. The asymmetry in the Hawkes model is an ongoing research topic, for example, consult El Euch et al. (2016). In our notation and setting, the asymmetric Hawkes model in El Euch et al. (2016) can be regarded as the Hawkes model of Eqs (1) and (2) with = ; = + ( 1) 12 21 22 11 21 volatility Mean signature plot Mean signature plot 0.35 φ = −0.5 κ = 0.5 0.30 φ = −0.3 0.10 κ = 1.0 φ = 0.0 0.25 κ = 2.0 φ = 0.1 0.08 φ = 0.3 0.20 0.06 0.15 0.04 0.10 0.02 0.05 0 0 60 120 180 240 300 60 120 180 240 300 Second Second Figure 5: Mean signature plot with xed = 0:5 with various (left) and xed = 0:3 and various (right) where is a newly introduced parameter. It is believed that there are many possible ways to incorporate asymmetry into the Hawkes model. In general, estimating in the di usion models is not trivial (Ait-Sahalia et al., 2013). One method for estimating is to use the method of moment as in Lee (2016). Let [X; Y ] denote the quadratic covariation process between the processes X and Y , i.e., Z Z t t [X; Y ] = X Y X dY Y dX t t t s s s s 0 0 = X Y + lim (X X )(Y Y ) 0 0 i+1 i i+1 i jj jj!0 for a sequence of random partitions with a limit in probability. The third moment variation of the return [R ; R] introduced by Choe and Lee (2014b) is a useful quantity to measure the skewness of the return distribution. In addition, the tractability of the di usion process enable us to easily derive the following formula. Proposition 5. Under the stationarity condition of the variance process with time 0, the following moment condition can be derived E[[R ; R] ] = K where 1 2 2 t 1 2 t 2 1 K = t 1 + e + t t + 1 e 2 1 3 2 3 S 2 0 2 2 1 2 2 2 2 2 ( )t + ( + )t 1 2 2 1 1 2 2 ( + ) 2 2 1 2 2 2 + + 1 2 1 2 1 2 1t 2t (1+2)t (e 1) (e 1) + (e 1) : 1 2 1 2 Proof. See F. If the drift part in the price process is zero, i.e., = 0, then the expectation of the third moment variation is simply 2 2 2 t E[[R ; R] ] = t 1 + e t t 2 3 2 3 S S 0 2 0 where the approximation is for a suciently large t. \ 2 [R ;R] i \ 2 Example 1. By Proposition 5, ! as the sample size increases where [R ; R] denotes N K the realized nite sum approximation of the third moment variation. The convergence of the estimates of in Figure 6 were plotted in a simulation study with parameter settings = 1:15; = 0:45; = 4 \ 2 2:8 10 ; = 0:85; = 0:0375; = 0:5. The sample mean of [R ; R] =K converges to . In the simulation result, the sample mean is 0:4935 with the standard error of 0:0757. However, it should be noted that the number of samples should be sucient for the convergence. If the number of samples is not sucient, it is better to use the approximate likelihood method or the simulated likelihood estimate discussed in 3.4. 13 1 −1 −2 −3 200 400 600 800 1000 sample size Figure 6: Convergence of the estimates of Figure 7: Hawkes model and di usion analogy 3.3 Comparison Both the Hawkes and the di usion models well describe the microstructure of price dynamics such as trade clustering or microstructure noise. The Hawkes model directly describes the tick-by-tick structure of the asset price and data is applied to the model without further assumptions or data corrections. The model's closed-form formula of the log-likelihood function and quite reliable numerical algorithms to nd the maximum the applicable. On the other hand, the di usion approach naturally extends the methodology traditionally used to describe asset price movements. Note that the di usion model in our paper is not a rigorous mathematical transform of the Hawkes model. We use the derivation to provide an intuition not a mathematical proof. Thus, one can argue about the legitimacy of the model, for example, the introduction of which we regarded as a constant. Nevertheless, the model inherits the advantages of typical di usion models. Based on the It o calculus and PDE approach, the derivations of useful formula such as moment conditions and distributional prop- erty are simpler than the Poisson based Hawkes models. Since the di usion model has been extensively studied for a long time, it is expected that there will be a more convenient aspect to apply the existing theory or extend the model. Meanwhile, the maximum likelihood estimation for the di usion model is generally more complicated because the closed-form formula for the density function is not available in many cases. In the absence of the closed-form likelihood function, the expansion based likelihood function approach (A t-Sahalia et al., 2008), simulation based method (Brandt and Santa-Clara, 2002) or the generalized method of moment (Garcia et al., 2011; Bollerslev et al., 2011) are used to estimate the parameters. 3.4 Simulated likelihood estimation Because the exact likelihood formula of the di usion process in this paper is barely available, the estima- tion is based on the simulation method proposed by Brandt and Santa-Clara (2002). Brie y explaining the method, the interval between two observed points, t and t , are divided into subintervals with a i i+1 length N . The M number of paths are simulated from t up to N 1 subintervals using the discretized version of the di usion model. The mean of the transition probability functions from the last values ρ Table 2: minimal tick percentage (%) - mid price BAC CVX GE IBM JPM KO MCD T VZ XOM 2007 79.91 61.88 83.29 55.70 78.16 75.41 73.71 79.15 84.57 70.14 2008 87.83 59.01 79.93 43.65 59.68 67.48 58.07 68.58 69.48 68.57 2009 88.19 70.36 93.84 56.98 72.42 80.10 84.51 82.14 82.06 86.79 2010 79.78 87.83 98.74 77.48 95.77 94.61 83.88 82.40 82.71 86.98 2011 99.53 72.16 99.21 52.73 96.96 89.44 86.92 98.23 89.07 90.35 Table 3: Minimal tick percentage (%) - transacted price BAC CVX GE IBM JPM KO MCD T VZ XOM 2007 92.70 69.86 97.26 64.40 90.34 89.81 87.23 94.74 93.11 79.78 2008 84.60 50.88 89.81 51.80 71.21 76.61 63.89 86.87 82.79 60.93 2009 98.89 72.34 98.39 59.98 88.03 89.15 80.42 97.44 93.67 82.47 2010 99.63 81.07 99.61 80.88 95.58 92.59 85.65 99.15 98.19 92.36 2011 99.81 62.08 99.68 57.04 96.76 91.76 80.58 99.01 97.12 84.48 of the simulated paths to the observed value at t , which is approximated by the normal distribution i+1 based on the discretization, becomes the maximum simulated likelihood. In Empirical studies, the data is reformulated to apply the di usion model because the original data is based on a tick structure. The large interval, i.e., t t , is set to one minute where a suciently large i+1 i number of events are observed for the approximation. Figure 7 presents the procedure, with every one minute, the observed price is the base point to construct a di usion process, which lies behind the tick structure. Within the interval, the paths of the discretized version of the di usion model are simulated with 60 subintervals. 4 Empirical study 4.1 Data For empirical studies, ultra high-frequency data of 10 stocks in the S&P 500 are used. As raw data in the rst place, we reorganize the data in the following way: The historical data consists of the best bid, ask quotes of the stocks, and their dynamics over trading time with various exchanges. The mid-price dynamics of the best bid and ask quotes of each stock reported in the New York Stock Exchange (NYSE) from 10:00 to 15:30 are selected to avoid the seasonal e ects observed in early or late in the market. In the original raw data, the time stamps have 1 second resolutions. If the prices changes are reported several times for one second, the price changes with equidistant intervals are redistributed over one second. The mid-price increments and decrements have a unit size of change that is the half of the minimal bid ask spread. If a price increment or decrement is larger than the minimal unit size, the change is considered to be the sum of the successive movements with the minimal size. In recent data, the percentage of the minimal change is very high in many symbols, as listed in Table 2. In addition, Table 3 lists the percentage of minimal change of transacted prices where similar patterns to the percentage of the mid-prices are observed. The symbols in the table represents: BAC - Bank of America Corp, CVX - Chevron, GE - General Electric Co., IBM - International Business Machines, JPM - JP Morgan Chase & Co, KO - The Coca-Cola Company, MCD - McDonald's Corp, T - AT&T Inc, VZ - Verizon Communications Inc, XOM - Exxon Mobil Corp 15 Table 4: Symmetric Hawkes estimation result, GE, January 2011 Date H.vol TSRV RRV s c 0103 0.0067 0.4661 1.3576 2.2596 0.0957 0.1289 0.1103 (0.0004) (0.0609) (0.0958) (0.1160) 0104 0.0082 0.4853 1.3941 2.5297 0.1139 0.1468 0.1344 (0.0005) (0.0494) (0.0889) (0.1104) 0105 0.0112 0.4741 1.1698 2.2281 0.1339 0.1825 0.1619 (0.0006) (0.0402) (0.0673) (0.0939) 0106 0.0091 0.5599 1.0112 2.1822 0.1265 0.1656 0.1391 (0.0005) (0.0471) (0.0675) (0.0958) 0107 0.0163 0.6973 0.5968 1.9959 0.1933 0.1998 0.1932 (0.0007) (0.0407) (0.0391) (0.0747) 0110 0.0132 0.4978 0.7434 1.8366 0.1520 0.1730 0.1553 (0.0006) (0.0360) (0.0465) (0.0780) 0111 0.0081 0.6959 0.7522 2.1448 0.1310 0.1414 0.1254 (0.0005) (0.0593) (0.0645) (0.1125) 0112 0.0098 0.4210 0.7322 1.8399 0.1181 0.1440 0.1328 (0.0005) (0.0406) (0.0559) (0.0965) 0113 0.0097 0.6512 0.3476 1.7471 0.1533 0.1275 0.1434 (0.0005) (0.0526) (0.0379) (0.1035) 0114 0.0092 0.6173 0.5303 1.8327 0.1390 0.1328 0.1344 (0.0005) (0.0533) (0.0498) (0.1097) 0118 0.0097 0.5122 0.4779 1.6361 0.1351 0.1235 0.1161 (0.0005) (0.0449 (0.0434) (0.0976) 0119 0.0168 0.4737 0.5202 1.5618 0.1774 0.1842 0.1772 (0.0007) (0.0331) (0.0351) (0.0772) 0120 0.0181 0.6944 0.5065 1.9431 0.2071 0.1913 0.1928 (0.0007) (0.0428) (0.0359) (0.0835) 0121 0.0316 0.5334 0.5986 1.7352 0.2358 0.2168 0.2299 (0.0010) (0.0262) (0.0282) (0.0556) 0124 0.0120 0.4355 0.4512 1.4181 0.1408 0.1249 0.1388 (0.0006) (0.0352) (0.0358) (0.0836) 0125 0.0209 0.6164 0.4259 1.6312 0.2144 0.2018 0.2096 (0.0008) (0.0321) (0.0275) (0.0640) 0126 0.0143 0.5053 0.5590 1.6587 0.1532 0.1146 0.1326 (0.0007) (0.0365) (0.0390) (0.0815) 0127 0.0147 0.5029 0.3925 1.4351 0.1687 0.1807 0.1578 (0.0007) (0.0371) (0.0317) (0.0868) 4.2 Dynamics of parameters and performance of volatility measure The parameters of the symmetric Hawkes process were estimated, as explained in C using the mid-price dynamics of the stocks quoted in NYSE. The estimations are employed on a daily basis because there are enough samples even in a day and the aim is to demonstrate the daily change in the parameters. Table 4 lists one of the results with GE for each day from January 3 to 27, 2011. The estimates of ; ; ; s c and their numerically computed standard errors in the parentheses are reported. In the estimation, the unit time, t = 1, is one second. The averaged daily estimates of ; ; ; for the di erent stocks over s c a month, January 2011, are also reported. In the column, `H.vol', the annualized daily volatility estimates computed by the estimates of the Hawkes parameters and using the formula in Proposition 3 are presented. In the column, `TSRV', the two scaled realized volatilities introduced by Zhang et al. (2005) are compared and in the column, `RRV', the volatility estimates proposed by Robert and Rosenbaum (2011) based on the uncertainty zones model are presented. The table shows that the Hawkes volatility, TSRV and RRV have similar values all over the reported time. Figure 8 plots the dynamics of the parameters of GE, 2011. The estimation results show evidence that the parameters of the Hawkes process, particularly , changes with time. The dynamics of with 16 Table 5: Averaged estimation result of symmetric Hawkes model, January 2011 Symbol s c GE 0.0141 0.5480 0.6766 1.8476 IBM 0.1489 1.0057 0.5037 2.0986 JPM 0.0672 0.6330 0.4767 1.5806 KO 0.0357 0.6669 0.3153 1.4814 MCD 0.0478 0.7201 0.4177 1.6641 T 0.0153 0.4593 0.5157 1.4506 VZ 0.0216 0.7887 0.4669 1.8206 XOM 0.0691 0.5280 0.3482 1.2808 time shows the typical movements of positively autocorrelated time series, which is strongly associated with the macro feature of the volatility movement, such as the GARCH e ect and stochastic volatility. In addition, a comparison of Figures 8a, 8b and 8c veri es that the dynamics of is related signi cantly to the dynamics of the volatility. When the parameter of a day is large, the volatility of the day is large and when the parameter of a day is small, the volatility of the day is small. In the gure, the movements of the other parameters ; and , do not appear to be meaningful s c compared to the movement of . The plots also show that the volatilities computed by the symmetric Hawkes modeling and TSRV show similar patterns over the observed time period. Figure 9 presents the parameter and volatility dynamics of GE, 2010. Similar to the previous case, the behaviors and TSRV are similar. The day of peaked volatilities in the gure is the May 6, 2010 Flash Crash where the equity prices fell rapidly. At the day of the Flash Crash, the two estimated volatilities had di erent values and TSRV was much larger than the Hawkes volatility. In addition, in Figure 10, the estimated Hawkes volatility, TSRV and RRV of T (left) and MCD (right) are compared. All three volatilities have similar forms of movements during the observed period. For T, the Hawkes volatility was close to TSRV (right) at the day of Flash Crash. The estimated Hawkes volatility of MCD at the Flash Crash was larger than the TSRV or RRV. Figure 11 plots the dynamics of the estimated parameters of GE in 2008, the starting year of the global nancial crisis. The dramatic changes in the Hawkes volatility, TSRV, ; and were observed in the beginning of the crisis around August 2008. The mutual excited parameter, , was rather stable. The volatility estimation results were compared using the Hawkes model and TSRV method in Table 6 with 10 symbols from 2007 to 2011. In each panel of the table, the mean of the Hawkes volatility and TSRV for given year and mean percentage error are presented. The volatility estimated by the symmetric Hawkes model is generally larger than the TSRV and the di erences between the two volatilities are around 15-25%. The reason for the discrepancy between the TSRV and the Hawkes volatility is unknown. Possible reasons include the intraday variation of the parameters as in the Flash Crash and the restrictions in the parameter condition for the symmetry. These two issues are examined in the following subsections. 4.3 Intraday volatility One of the interesting applications to modeling the daily price dynamics using the symmetric Hawkes process is that the intraday volatility can be estimated in almost every moment of the day. This is possible because every arrival time of price change, which are plentiful even during ten minutes, is used and the maximum likelihood estimation is so powerful that the parameters can be estimated with similar or less than ten minutes data. Figure 12 shows the dynamics of the intraday volatility of GE with randomly chosen days. The rst estimation of each day was performed using the rst ten minutes data of each day. In this example, it ranged from 10:00 a.m. to 10:10 a.m. The price movement histories were then updated in every ten minutes and the intraday volatil- ities were re-estimated using the updated data and already existing one. For the estimation, the reparametrization in Remark 2 were used and hence the annualized volatility was estimated directly by the maximum likelihood estimation with its numerically computed standard error. The solid lines in the gure represent the annualized volatilities estimated by the intraday data up to the time and the dotted lines represent the standard errors. In Figures 12a and 12b, the volatilities are generally large at the beginning of the day and tend to decreases, which is consistent with the seasonality e ect in that in the early markets, more trading activities are observed than the middle of the day. 17 GE, 2011 GE, 2011 0.8 0.8 0.4 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (a) volatility (b) TSRV GE, 2011 GE, 2011 1.2 0.1 0.6 0 0 3 6 9 12 3 6 9 12 Month Month (c) (d) GE, 2011 GE, 2011 1.2 0.6 0 0 3 6 9 12 3 6 9 12 Month Month (e) (f ) Figure 8: Symmetric Hawkes estimation result, GE, 2011 Table 6: Comparison of the volatility estimation by the Hawkes model and realized volatility BAC CVX GE IBM JPM KO MCD T VZ XOM H.vol 0.1725 0.2387 0.1406 0.1674 0.2170 0.1290 0.1494 0.1731 0.1615 0.2327 TSRV 0.1555 0.1834 0.1323 0.1380 0.1899 0.1118 0.1296 0.1615 0.1433 0.1730 MPE(%) 14.93 20.28 12.34 17.26 15.91 15.58 15.87 14.04 14.26 23.14 H.vol 0.7008 0.5129 0.3646 0.3939 0.6929 0.2583 0.3286 0.3751 0.3806 0.4592 TSRV 0.4880 0.3077 0.3172 0.2658 0.4868 0.2031 0.2440 0.3069 0.2930 0.2827 MPE(%) 28.14 33.17 13.57 28.41 26.61 19.60 21.34 16.09 19.97 32.53 H.vol 0.7397 0.2608 0.3439 0.2166 0.4701 0.1762 0.1883 0.2407 0.2029 0.2198 TSRV 0.5420 0.2029 0.3367 0.1664 0.3722 0.1509 0.1607 0.1934 0.1810 0.1773 MPE(%) 25.42 21.48 9.12 20.81 17.97 15.23 15.79 20.27 13.43 18.26 H.vol 0.2869 0.1758 0.1963 0.1395 0.2223 0.1138 0.1258 0.1461 0.1578 0.1603 TSRV 0.2234 0.1376 0.1952 0.1184 0.1985 0.1031 0.1040 0.1255 0.1241 0.1291 MPE(%) 21.20 21.11 11.47 17.07 13.46 13.51 18.12 15.96 17.73 18.86 H.vol 0.3426 0.2197 0.2389 0.1681 0.2554 0.1334 0.1297 0.1460 0.1543 0.1872 TSRV 0.2648 0.1719 0.1921 0.1334 0.2146 0.1091 0.1122 0.1227 0.1228 0.1553 MPE(%) 13.79 20.46 18.78 18.49 17.06 19.24 15.73 17.47 20.38 17.38 α μ volatility TSRV β GE, 2010 GE, 2010 1.6 1.6 1.2 1.2 0.8 0.8 0.4 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (a) volatility (b) TSRV GE, 2010 GE, 2010 0.1 0.5 0 0 3 6 9 12 3 6 9 12 Month Month (c) (d) GE, 2010 GE, 2010 1.8 1.2 0.6 0 0 3 6 9 12 3 6 9 12 Month Month (e) (f ) Figure 9: Symmetric Hawkes estimation result, GE, 2010 α μ volatility TSRV β T, 2010 MCD, 2010 0.8 0.8 0.4 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (a) T, volatility (b) MCD, volatility T, 2010 MCD, 2010 0.8 0.8 0.4 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (c) T, TSRV (d) MCD, TSRV T,2010 MCD,2010 0.8 0.8 0.4 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (e) T, RRV (f ) MCD, RRV Figure 10: Volatility comparisons with symmetric Hawkes estimation results, T (left) and MCD (right), RRV TSRV volatility RRV TSRV volatility GE, 2008 GE, 2008 1.6 1.6 1.2 1.2 0.8 0.8 0.4 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (a) volatility (b) TSRV GE, 2008 GE, 2008 0.4 1.6 0.3 1.2 0.2 0.8 0.1 0.4 0 0 3 6 9 12 3 6 9 12 Month Month (c) (d) GE, 2008 GE, 2008 0.6 3 0.4 2 0.2 1 0 0 3 6 9 12 3 6 9 12 Month Month (e) (f ) Figure 11: Symmetric Hawkes estimation result, GE, 2008 α μ volatility TSRV β GE, 2011/08/05 GE, 2011/10/13 1.5 0.8 0.6 0.4 0.5 0.2 0 0 10:00 12:00 14:00 15:30 10:00 12:00 14:00 15:30 (a) (b) GE, 2010/05/06 VZ, 2010/05/06 1.0 1.0 0.8 0.8 0.6 0.6 0.4 0.4 0.2 0.2 0 0 10:00 12:00 14:00 15:30 10:00 12:00 14:00 15:30 (c) (d) Figure 12: Estimated cumulative intraday volatility (annualized) with every ten minutes update Figure 12c shows the data for the 2010/05/06 Flash Crash and a dramatic increase was observed in the late part of the day. Similar behavior is presented in Figure 12d, which is for the intraday volatility of VZ in the 2010/05/06 Flash Crash. The real time volatility measurement technique will be very useful in intraday risk managements, because investors can respond to sudden market changes more e ectively, if they can compute the exact volatility variation. 4.4 Fully characterized Hawkes process The maximum likelihood estimation of the Hawkes process with full characterization of the parameters as explained in Subsection 2.2 was performed. As in the previous subsection, the estimations were employed on a daily basis. In Figure 13, the dynamics of all parameters of the fully characterized Hawkes model of GE in 2011 are plotted. The dynamics of the parameters ; are close to each other in the mean, as illustrated in Figure 13a 1 2 which suggest that = in the long run sense. Similarly, each pair of parameters of ( ; ), 1 2 11 22 ( ; ), ( ; ), and ( ; ) are close to each other in the mean. The dynamics of the parameters 12 21 11 22 12 21 and uctuate more than and over time. The sample means of and of GE in 2011 12 21 11 22 12 21 are quite close to and as reported in panel A of Table 7. The row `std.' in the table is the sample 11 22 standard deviation of the time series of each parameter over the time period. On the other hand, the parameters are not always close to among others in the mean. In the ij panel B of the table which presents the estimates of the parameters of XOM in 2008, the estimates of and are close to each other and similarly, the estimates of and are close to each other in the 22 12 21 mean, but the estimates of and are signi cantly di erent in the mean. Similarly, the di erence 11 12 in the estimates of and are signi cant. In this case, the self-excited e ects are less persistent 21 22 than the mutually excited e ects. 4.5 Di usion parameter The parameters of the di usion model introduced in subsection 3.1 were estimated using the simulated likelihood estimation explained in subsection 3.4. The results of GE, January 2011 are presented in Table 8 and are similar to the results with the Hawkes model in Table 4. In addition, the estimates of the di usion model with are presented in Table 9. The estimates in each model show the similar patterns over the period. The di usion estimation has its own pros and cons. In our setting, because the observed values of the price over one minute intervals are only used, which is in contrast to the Hawkes modeling where all times volatility volatility volatility volatility GE, 2011 1.2 0.08 0.04 3 6 9 12 Month (a) and 1 2 GE, 2011 1.5 0.5 3 6 9 12 Month (b) and 11 22 GE, 2011 1.5 0.5 3 6 9 12 Month (c) and 12 21 JPM, 2011 3 6 9 12 Month (d) and 11 22 GE, 2011 3 6 9 12 Month (e) and 12 21 Figure 13: Estimation result with the fully characterized Hawkes, GE, 2011 α α μ β β c s c s Table 7: Estimation result of fully characterized self and mutually excited Hawkes process, GE, 2011 in panel A and XOM, 2008 in panel B 1 2 11 22 12 21 11 22 12 21 mean 0.0198 0.0199 0.5196 0.5228 0.3235 0.3165 1.4145 1.4128 1.5574 1.5378 std. 0.0124 0.0128 0.1108 0.1241 0.1902 0.1875 0.2097 0.2463 0.5856 0.5402 mean 0.1886 0.1727 1.0594 0.9904 0.1369 0.1288 1.7648 1.6916 0.8972 0.7329 std. 0.1187 0.1196 0.3274 0.3348 0.0850 0.0786 0.4478 0.4480 0.6104 0.5197 Table 8: Di usion model estimation result, GE, January 2011 Date m a a b volatility s c 0103 0.0137 0.1225 1.9872 2.5512 0.1392 0104 0.0091 0.0805 2.9568 3.3122 0.1546 0105 0.0089 0.7037 2.7975 3.5694 0.3661 0106 0.0082 1.0110 2.1844 3.3689 0.2511 0107 0.0196 0.8574 1.3001 2.5649 0.2557 0110 0.0140 0.7615 0.7941 2.3103 0.1739 0111 0.0069 0.5888 2.1067 2.9101 0.1710 0112 0.0046 0.3888 1.5751 2.4639 0.1247 0113 0.0105 0.6536 0.7633 2.3407 0.1314 0114 0.0043 0.9421 1.2195 2.2284 0.2833 0118 0.0093 0.7737 0.7164 1.7152 0.2325 0119 0.0174 0.5661 1.1345 2.3551 0.1714 0120 0.0173 0.7892 1.0437 2.6670 0.1844 0121 0.0342 0.5358 1.1074 2.3377 0.2215 0124 0.0142 0.4669 0.8490 2.0244 0.1349 0125 0.0208 0.7222 0.9767 2.1042 0.2308 0126 0.0135 0.6058 1.0571 2.4183 0.1387 0127 0.0157 0.6458 0.8145 1.8174 0.2043 of price changes are used, the di usion estimator is less ecient than the Hawkes estimator. In addition, by the nature of the simulated likelihood estimation, it takes longer time to compute the likelihoods and the computed results are not deterministic but depend on the random numbers generated by computers. On the other hand, when the observing times of a price process are limited, i.e., the prices are only available at each one minute interval, the di usion model and its estimation are a feasible alternate choice to examine the nature of the price movements in high-frequency. 5 Conclusion This paper examined the empirical performance of the symmetric Hawkes process which is a simple model to consider for both clustering property and market microstructure noise in volatility estimation using the stock prices in the S&P 500. The daily dynamics of the Hawkes parameters, the comparison between the Hawkes volatility and the realized volatility and the intraday volatility estimation procedure are discussed. The di usion analogy of the symmetric Hawkes model was also proposed to provide the analytical simplicity for computing the distributional properties. The di usion model also incorporates the clustering e ect, market microstrucutre noise, in addition to asymmetric property. The volatility could be estimated over a relatively short time interval with the Hawkes model and the intraday variations of volatility was demonstrated. A comparison between the Hawkes volatility and TSRV showed the di erence around 15-25%. The parameter restriction, asymmetry and parameter variations might be the cause of the discrepancy but more work will be needed to understand the exact reason. The estimation results of the di usion model were provided where similar patterns to the Hawkes model parameters were observed. 24 Table 9: Di usion model estimation result with , GE, January 2011 Date m a a b s c 0103 0.0140 0.0191 1.9107 2.4863 0.2010 0104 0.0094 0.4947 2.3869 3.2896 0.1106 0105 0.0083 0.4363 3.0567 3.5328 -0.0338 0106 0.0097 0.7811 2.3254 3.4391 -0.0952 0107 0.0196 0.8110 1.5196 2.6318 -0.3633 0110 0.0067 0.5657 1.7872 2.3712 0.0397 0111 0.0072 0.7457 1.7410 2.5917 0.0229 0112 0.0116 0.4343 1.4936 2.7054 -0.1048 0113 0.0114 0.7440 0.8413 2.2652 -0.0709 0114 0.0051 0.8205 1.2866 2.1705 -0.0087 0118 0.0045 0.1585 1.6107 1.7917 0.1269 0119 0.0185 0.4916 1.1835 2.3827 -0.2130 0120 0.0153 0.5343 0.9029 1.9994 0.1086 0121 0.0164 0.6817 1.1065 1.9896 -0.1220 0124 0.0148 0.5293 1.0546 1.9986 -0.0055 0125 0.0252 0.6536 0.9445 2.0343 0.1237 0126 0.0156 0.5248 1.0506 2.3678 0.0473 0127 0.0167 0.5304 0.8537 1.6645 0.0800 References A t-Sahalia, Y., Cacho-Diaz, J., and Laeven, R. J. (2010). Modeling nancial contagion using mutually exciting jump processes. Technical report, National Bureau of Economic Research. A t-Sahalia, Y. et al. (2008). Closed-form likelihood expansions for multivariate di usions. The Annals of Statistics, 36(2):906{937. Ait-Sahalia, Y., Fan, J., and Li, Y. (2013). The leverage e ect puzzle: Disentangling sources of bias at high frequency. Journal of Financial Economics, 109:224{249. A t-Sahalia, Y., Mykland, P. A., and Zhang, L. (2005). How often to sample a continuous-time process in the presence of market microstructure noise. Review of Financial studies, 18:351{416. A t-Sahalia, Y., Mykland, P. A., and Zhang, L. (2011). Ultra high frequency volatility estimation with dependent microstructure noise. Journal of Econometrics, 160:160{175. Andersen, T. G., Bollerslev, T., Diebold, F. X., and Labys, P. (2003). Modeling and forecasting realized volatility. Econometrica, 71:579{625. Bacry, E., Dayri, K., and Muzy, J.-F. (2012). Non-parametric kernel estimation for symmetric Hawkes processes. application to high frequency nancial data. The European Physical Journal B-Condensed Matter and Complex Systems, 85:1{12. Bacry, E., Delattre, S., Ho mann, M., and Muzy, J.-F. (2013). Modelling microstructure noise with mutually exciting point processes. Quantitative Finance, 13:65{77. Bacry, E. and Muzy, J.-F. (2014). Hawkes model for price and trades high-frequency dynamics. Quan- titative Finance, 14:1147{1166. Barndor -Nielsen, O. E. and Shephard, N. (2002a). Econometric analysis of realized volatility and its use in estimating stochastic volatility models. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 64:253{280. Barndor -Nielsen, O. E. and Shephard, N. (2002b). Estimating quadratic variation using realized vari- ance. Journal of Applied Econometrics, 17:457{477. Bauwens, L. and Hautsch, N. (2009). Modelling nancial high frequency data using point processes. In Handbook of Financial Time Series, pages 953{979. Springer Berlin Heidelberg. 25 Bollerslev, T., Gibson, M., and Zhou, H. (2011). Dynamic estimation of volatility risk premia and investor risk aversion from option-implied and realized volatilities. Journal of econometrics, 160(1):235{245. Bowsher, C. G. (2007). Modelling security market events in continuous time: Intensity based, multivariate point process models. Journal of Econometrics, 141:876{912. Brandt, M. W. and Santa-Clara, P. (2002). Simulated likelihood estimation of di usions with an appli- cation to exchange rate dynamics in incomplete markets. Journal of nancial economics, 63:161{210. Br emaud, P. (1981). Point Processes and Queues : Martingale Dynamics. Springer. Broyden, C. G. (1970). The convergence of a class of double-rank minimization algorithms 1. general considerations. IMA Journal of Applied Mathematics, 6:76{90. Choe, G. H. and Lee, K. (2014a). Conditional correlation in asset return and garch intensity model. AStA Advances in Statistical Analysis, 98:197{224. Choe, G. H. and Lee, K. (2014b). High moment variations and their application. Journal of Futures Markets, 34:1040{1061. Da Fonseca, J. and Zaatour, R. (2014a). Clustering and mean reversion in a Hawkes microstructure model. Journal of Futures Markets, 35:813{838. Da Fonseca, J. and Zaatour, R. (2014b). Hawkes process: Fast calibration, application to trade clustering, and di usive limit. Journal of Futures Markets, 34:548{579. Daley, D. J. and Vere-Jones, D. (2003). An introduction to the theory of point processes, volume 1. Springer. Dassios, A. and Zhao, H. (2012). Ruin by dynamic contagion claims. Insurance: Mathematics and Economics, 51:93{106. El Euch, O., Masaaki, F., and Mathieu, R. (2016). The microstructural foundations of leverage e ect and rough volatility. arXiv preprint arXiv:1609.05177. Embrechts, P., Liniger, T., Lin, L., et al. (2011). Multivariate hawkes processes: an application to nancial data. Journal of Applied Probability, 48:367{378. Errais, E., Giesecke, K., and Goldberg, L. R. (2010). Ane point processes and portfolio credit risk. SIAM Journal on Financial Mathematics, 1:642{665. Garcia, R., Lewis, M.-A., Pastorello, S., and Renault, E. (2011). Estimation of objective and risk-neutral distributions based on moments of integrated volatility. Journal of Econometrics, 160(1):22{32. Hansen, P. R. and Lunde, A. (2006). Realized variance and market microstructure noise. Journal of Business and Economic Statistics, 24:127{161. Hawkes, A. G. (1971a). Point spectra of some mutually exciting point processes. Journal of the Royal Statistical Society. Series B (Methodological), pages 438{443. Hawkes, A. G. (1971b). Spectra of some self-exciting and mutually exciting point processes. Biometrika, 58:83{90. Hawkes, A. G. and Oakes, D. (1974). A cluster process representation of a self-exciting process. Journal of Applied Probability, pages 493{503. Henningsen, A. and Toomet, O. (2011). maxlik: A package for maximum likelihood estimation in R. Computational Statistics, 26:443{458. Heston, S. (1993). A closed-form solution for options with stochastic volatility with applications to bond and currency options. Review of Financial Studies, 6:327{343. Hewlett, P. (2006). Clustering of order arrivals, price impact and trade path optimisation. In Workshop on Financial Modeling with Jump processes, Ecole Polytechnique. 26 Jaisson, T., Rosenbaum, M., et al. (2015). Limit theorems for nearly unstable Hawkes processes. The Annals of Applied Probability, 25:600{631. Large, J. (2007). Measuring the resiliency of an electronic limit order book. Journal of Financial Markets, 10:1{25. Lee, K. (2016). Probabilistic and statistical properties of moment variations and their use in inference and estimation based on high frequency return data. Studies in Nonlinear Dynamics & Econometrics, 20:19{36. Li, D.-H. and Fukushima, M. (2001). On the global convergence of the bfgs method for nonconvex unconstrained optimization problems. SIAM Journal on Optimization, 11:1054{1064. Nash, J. C. et al. (2014). On best practice optimization methods in R. Journal of Statistical Software, Ogata, Y. (1978). The asymptotic behaviour of maximum likelihood estimators for stationary point processes. Annals of the Institute of Statistical Mathematics, 30:243{261. Robert, C. Y. and Rosenbaum, M. (2011). A new approach for the dynamics of ultra-high-frequency data: The model with uncertainty zones. Journal of Financial Econometrics, 9:344{366. Zhang, L., Mykland, P. A., and A t-Sahalia, Y. (2005). A tale of two time scales. Journal of the American Statistical Association, 100:1394{1411. Zheng, B., Roue , F., and Abergel, F. (2014). Modelling bid and ask prices using constrained Hawkes processes: Ergodicity and scaling limit. SIAM Journal on Financial Mathematics, 5:99{136. A Expected intensity of fully characterized Hawkes model Consider the conditional expectation of the intensity processes: ` (tjs) = E[ (t)jF ]; ` (tjs) = E[ (t)jF ]: i i s ij ij s Then, for each i, Z Z s t (t s) (s u) (t u) ii ii ii ` (tjs) = E e e dN (u) + e dN (u) F ii ii i ii i s 1 s (t s) (t u) ii ii = (s)e + E e dN (u) F ii ii i s Z Z t t (t s) (t u) (t u) ii ii ii = (s)e + E e (dN (u) (u)du) + e (u)du F ii ii i i ii i s s s (t s) (t u) ii ii = (s)e + e ` (ujs)du ii ii i and by di erentiating both sides with respect to t, d` (tjs) ii (t s) (t u) ii ii = (s) e + ` (tjs) e ` (ujs)du ii ii ii i ii ii i dt = ` (tjs) ` (tjs) ii i ii ii = + ( )` (tjs) + ` (tjs): ii i ii ii ii ii ij In addition, with the similar method, for i 6= j, d` (tjs) ij = ( + ` (tjs) ` (tjs)) ` (tjs): ij j ji jj ij ij dt 27 The di erential equation system is represented by the matrix form 2 3 2 32 3 2 3 ` (tjs) 0 0 ` (tjs) 11 11 11 11 11 1 6 7 6 76 7 6 7 ` (tjs) 0 ` (tjs) 12 12 12 12 12 2 6 7 6 76 7 6 7 = + : 4 5 4 54 5 4 5 ` (tjs) 0 ` (tjs) 21 21 21 21 21 1 ` (tjs) 0 0 ` (tjs) 22 22 22 22 22 2 When the eigenvalues of the matrix are negative, the particular solution of the system becomes the long-run expectations of the intensities and are given by 2 3 2 3 ` (tjs) f( ) + g 11 11 21 22 22 12 1 12 22 2 6 7 6 7 ` (tjs) 1 f( ) + g 12 12 22 11 11 21 2 21 11 1 6 7 6 7 4 5 4 5 ` (tjs) f( ) + g 21 21 11 22 22 12 1 11 22 2 ` (tjs) f( ) + g 22 22 12 11 11 21 2 21 11 1 as t ! 1, where H = ( ) ( + ): 11 12 21 22 22 11 22 12 21 12 21 22 12 21 22 The above formulas can be used as the presumed initial values of the intensity processes in the simulations or estimation procedures. With full characterization of the parameters ; ; , the system is four ij ij ij dimensional and the solution is rather complicated. B Simulation method If the decaying parameters are di erent from each other, the system of the self and mutually excited ij Hawkes and intensity processes (N ; N ; ; ) are not Markov. As shown in Eq. (1), (t) depends 1 2 1 2 1 on both (s) and (s), for s < t, and similarly, (t) depends on both (s) and (s). On the 11 12 2 21 22 other hand, the whole system of the processes (N ; N ; ; ; ; ) are Markov and to generate the 1 2 11 12 21 22 future paths, it is only important to know the current values of (N ; N ; ; ; ; ) not the entire 1 2 11 12 21 22 past histories of the processes. Therefore, for the simulation of the Hawkes process, it is important to compute the distributions of the arrival times determined by each component of the intensities, and ij Suppose that, over a time interval [s; t), there is no jump by N and N ; then the intensities are 1 2 deterministic and exponentially decaying function is (t s) ij (t) = (s)e : ij ij Note that N (t) N (s) can be represented by the sum of three jump components N ; N , and N i i i0js iijs ijjs independent upon F with the corresponding intensities , , and , respectively. Let be the rst s i ij ij ijjs interarrival time of N with intensity after s. The probability distribution of is then represented ij ij ijjs by ij 1 e Pf > ug = exp (s) : ijjs ij ij Thus, 1 log U ij log 1 + ijjs (s) ij ij where U is a uniformly distributed random variable over [0,1]. In addition, let denote a random i0 variable that follows a Poisson distribution with intensity . Then minf ; ; g becomes the next i 10 20 ijjs jump arrival time after s. After a jump occur, the counting processes are updated accordingly, the intensities are updated, as in Eqs. (1) and (2), and the above procedure is applied repeatedly. C Likelihood function Let t be the k-th jump arrival time of N and be the interarrival time between k and (k + 1)-th k 1 1jk jumps. Then the conditional cumulative distributions of at time t , i.e., with given (t ), is k 1 k 1jk t +u F (uj (t )) = 1 exp (s)ds : 1 k 1 1jk 28 Therefore, the conditional density functions is t +u f (uj (t )) = (t ) exp (s)ds : 1 k 1 k 1 1jk Similarly, let t be the m-th jump arrival time of N and be the interarrival time between the m m 2 2jm and (m + 1)-th jumps. The conditional density function of at time t is then 2jm t +u f (uj (t )) = (t ) exp (s)ds : 2 m 2 m 2 2jm Now consider the interval [0; T ] over which the jumps are observed. The log-likelihood of the realized jump arrivals up to time T is represented by the sum of log-likelihood of all realized arrivals of N and N . That is Z Z t t k+1 m+1 X X L(; T ) = log (t ) exp (s)ds + log (t ) exp (s)ds 1 k 1 2 m 2 t t k m k m Z Z Z T T T = log (; t)dN (t) + log (; t)dN (t) ( (; t) + (; t))dt 1 1 2 2 1 2 0 0 0 where = f ; ; g denotes the parameter vector. The maximum likelihood estimator, , is the ij ij ij estimator which maximize L under the observations of realized jump arrivals of N and N . 1 2 De ne a matrix I () with each element " # 1 @ @ 1 @ @ 1 1 2 2 I () = E + dt : ij @ @ @ @ 0 1 i j 2 i j The maximum likelihood estimator converges to the true parameter value asymptotically normally in distribution with an asymptotic variance-covariance matrix I ( ), see Ogata (1978). For the maximum likelihood estimation in the statistical package R, consult Henningsen and Toomet (2011). D Proof of the variance formula in Proposition 3 In this section, the variance formula is derived under the symmetric Hawkes process assumption of the price process. When the price follows Eq. (4) with symmetric Hawkes process, the variance of the return is represented by Var(N (t) N (t) (N (0) N (0)): 1 2 1 2 S (0) To compute the variance of the return over time interval [0; t], the following results are needed. The intensities and are assumed to be in the stationary state at time 0. Under the assumption, the 1 2 variance of the price process is derived using the stochastic integration theory. The quadratic variation of X is de ned by [X ] = X 2 X dX t s s and the quadratic covariation of X and Y is de ned as Z Z t t [X; Y ] = X Y X dY Y dX : t t t s s s s 0 0 When the processes are quadratic pure jump processes, i.e., the quadratic (co)variation of the continuous part is zero, X X 2 2 [X ] = X + (X ) ; [X; Y ] = X Y + (X Y ): t 0 s t 0 0 s s 0<st 0<st Without a loss of generality, it is assumed that N (0) = N (0) = 0 in this proof. The next lemma is 1 2 stated without proof. 29 Lemma 6. Under the stationarity condition of the intensities at time 0, (a) E[ (t)] = E[ (t)] = (0) = (0) 1 2 1 2 (b) E[N (t)] = E[N (t)] = (0)t 1 2 1 (c) E [[N ] ] = E [[N ] ] = E [N (t)] = (0)t 1 2 1 1 t t 2 2 2 (d) E [[ ] ] = E [[ ] ] = (0) + ( + ) (0)t 1 2 1 1 s c t t (e) E [[ ; ] ] = (0) + 2 (0)t 1 2 s c 1 (f ) E [[N ; ] ] = E [[N ; ] ] = E[N (t)] = (0)t 1 1 2 2 s 1 s 1 t t (g) E [[N ; ] ] = E [[N ; ] ] = E[N (t)] = (0)t 1 2 2 1 c 1 c 1 t t Recall that s c M = : c s Lemma 7. Under the stationarity condition of the intensities at time 0, 2 2 2 E[ (t)] 1 1 + + 2 1 2 t 2 t 1 s c 1 2 = c e + c e (0)M 1 2 1 E[ (t) (t)] 1 1 2( + ) 1 2 2 s c for some constant c and c . 1 2 Proof. Note that E[ (t)] = E [[ ] ] + 2E (u)d (u) 1 1 1 1 t 2 2 2 2 = (0) + ( + ) (0)t + 2E (u) + ( ) (u) + (u) (u) du 1 1 s c 1 2 1 s c 1 2 2 2 2 = (0) + ( + + 2 ) (0)t + 2 ( )E [ (u)] + E [ (u) (u)] du 1 s s c s 1 2 1 s c 1 and Z Z t t E[ (t) (t)] = E [[ ; ] ] + E (u)d (u) + E (u)d (u) 1 2 1 2 1 2 2 1 0 0 2 2 = (0) + 2 (0)t + 2E (u) + (u) + ( ) (u) (u) du s c 1 1 c s 1 2 1 1 2 2 = (0) + 2( + ) (0)t + 2 E[ (u)] + ( )E[ (u) (u)] du: s c 1 c s 1 2 1 1 Therefore, a system of equations can be derived: 2 3 dE[ (t)] 2 2 2 6 7 E[ (t)] + + 2 s c 1 s c dt = 2 + (0) : 4 5 1 dE[ (t) (t)] E[ (t) (t)] 2( + ) 1 2 c s 1 2 s c dt The particular solution of the system is 2 3 2 2 2 + ( + ) + ( ) 2 + c s s c s 2 2 6 7 2 2 1 1 + + 2 ( ) 1 s c s 6 c 7 (0)M = (0) 1 1 3 2 4 5 2( + ) + 2 ( ) + 2 + 2 2 s c 2 s c s c s 2 2 ( ) where the inverse matrix of M is represented by " # 2 2 2 2 1 ( ) ( ) s c 1 s s c c M = = : c s 2 2 2 2 1 2 ( ) ( ) s s c c 30 In addition, the general solution is 2 2 2 E[ (t)] 1 1 + + 2 2 t 2 t 1 1 1 2 s c = c e + c e (0)M 1 2 1 E[ (t) (t)] 1 1 2( + ) 1 2 s c and with the initial condition of (8), 2 2 (0)( ) (0)( + ) 1 s c 1 s c c = ; c = : 1 2 4 4 1 2 Lemma 8. Under the stationary state condition of the intensities at time 0, we have c c E[ (t)N (t)] 1 1 1 1 1 1 1 2 t t 2 t 2 t 1 2 1 2 = d e + d e + e + e 1 2 E[ (t)N (t)] 1 1 1 1 1 2 1 2 2 2 1 + 1 1 1 2 s c (0) M t + M (M ) : 1 2 c 2 s c Proof. Note that Z Z t t E[ (t)N (t)] = E [[ ; N ] ] + E (u )dN (u) + E N (u)d (u) 1 1 1 1 1 1 1 1 0 0 = (0)t + E (u)du s 1 + f E[N (u)] + ( )E[ (u)N (u)] + E[ (u)N (u)]g du 1 s 1 1 c 2 1 Z Z t t = (0)t + E[ (u)]du + (0)udu s 1 1 0 0 + f( )E[ (u)N (u)] + E[ (u)N (u)]g du s 1 1 c 1 2 and Z Z t t E[ (t)N (t)] = E [[ ; N ] ] + E (u )dN (u) + E N (u)d (u) 1 2 1 2 1 2 2 1 0 0 = (0)t + E (u) (u)du c 1 1 2 + f E[N (u)] + ( )E[ (u)N (u)] + E[ (u)N (u)]g du 2 s 1 2 c 2 2 Z Z t t = (0)t + E[ (u) (u)]du + (0)udu c 1 1 2 1 0 0 + f( )E[ (u)N (u)] + E[ (u)N (u)]g du: s 1 2 c 1 1 In the above, E[ (u)N (u)] is replaced with E[ (u)N (u)] due to the symmetry of the model and 2 1 1 2 similarly, E[ (u)N (u)] is replaced with E[ (u)N (u)]. Therefore, following system of equations can 2 2 1 1 be derived: 2 3 dE[ (t)N (t)] 1 1 E[ (t)N (t)] (0) + E[ (t)] + (0)t 6 7 1 1 s 1 1 dt = M + 4 5 dE[ (t)N (t)] E[ (t)N (t)] (0) + E[ (t) (t)] + (0)t 1 2 1 2 c 1 1 2 1 dt E[ (t)N (t)] 1 1 1 1 2 t 2 t 1 2 = M + c (0) e + c (0) e 1 1 2 1 E[ (t)N (t)] 1 1 1 2 2 2 + + 2 1 s c + (0) t + M 2( + ) c 2 s c 31 where the previous lemma is used. The particular solution is A B 1 1 1 1 2 t 2 t 1 2 t + + k e + k e 1 2 A B 1 1 2 2 where A 1 1 1 = (0) M A 1 and 2 2 B + 1 s 1 1 2 s c = (0) M (M ) B 2 2 c 2 s c and k = c = ; k = c = . The general solution is 1 1 1 2 2 2 c c E[ (t)N (t)] 1 1 1 1 A B 1 2 1 1 t t 2 t 2 t 1 1 1 2 1 2 = d e + d e + e + e + t + 1 2 E[ (t)N (t)] 1 1 1 1 A B 1 2 2 2 1 2 and with the initial condition, (0)( ) (0)( + ) 1 s c 1 s c d = ; d = : 1 2 2 2 2 2 1 2 Proposition 9. Under the stationary state condition of the intensities at time 0, 2d 2d E[N (t)] 1 1 1 2 1 t t 1 2 = (e 1) + (e 1) E[N (t)N (t)] 1 1 1 2 1 2 c c 1 1 2 1 2 t 2 t 1 2 + (e 1) + (e 1) 2 2 1 1 1 2 2 2 1 + 1 1 2 1 s 1 2 s c (0) M t + 2M (M ) t : 1 2 0 c s c In addition, 2 2 t e 1 E[(N (t) N (t)) ] = 2 (0) t 1 : 1 2 1 2 2 1 1 Proof. The following is used: E[N (t)] = E [[N ] ] + 2E N (u )dN (u) 1 1 1 1 t = (0)t + 2E (u)N (u)du 1 1 1 and Z Z t t E[N (t)N (t)] = E [[N ; N ] ] + E N (u )dN (u) + E N (u )dN (u) 1 2 1 2 1 2 2 1 0 0 = 2E (u)N (u)du 1 2 along with the previous lemmas. In the above equation, E [[N ; N ] ] = 0 as the probability of the 1 2 simultaneous jumps of N and N is zero. In addition, E[ (u)N (u)] = E[ (u)N (u)] under the 1 2 1 2 2 1 symmetry. Therefore, 2 2 E[(N (t) N (t)) ] = 2(E[N (t)] E[N (t)N (t)]) 1 2 1 2 8d 4c 2 1 1 1 t 2 t 2 1 1 = (e 1) (e 1) 2 (0) ( ) ( ) 1 t 1 s c s c 2 2 1 1 1 1 8d 4c 1 1 s c t 2 t 1 1 = (e 1) (e 1) + 2 (0) 1 t 1 1 2 2 4( ) ( ) s c s c t 2 t 1 1 = (0) (e 1) + (0) (e 1) + 2 (0) t 1 1 1 3 3 2 1 1 1 t 2 t 1 1 2 (0) e 1 e 1 2 2 = t 2( ) + ( ) : s c s c 1 1 32 E Proof of Proposition 4 Note that " # Z Z t t E n du + V dW u u u 0 0 " # " # Z Z Z Z 2 2 t t t t p p s s = E n du + 2E n du V dW + E V dW u u u u u u 0 0 0 0 " # Z Z Z Z t t t t = E n du + 2E n du V dW + E V du u u u u 0 0 0 0 " # Z Z Z t t t = E n du + 2E n du V dW + t: u u u 0 0 0 Let " # Z Z Z t t s x(t) := E n du = 2E n n du ds u s u 0 0 0 and Z Z t t y(t) :=E n du V dW u u 0 0 Z Z Z Z t s t s p p s s =E n du V dW + E V dW n ds u s u s s u 0 0 0 0 Z Z t s =E V dW n ds : u s 0 0 By assumming n = 0, Z Z Z Z Z s s s s s p p p p s s s s E n V dW = E n du V dW + E V dW V dW s u 1 u u u u u u u u 0 0 0 0 0 = y(s) + E[V ]du 1 u = y(s) + s: Thus, y (t) = y(t) + t and ( t 1 + e ) y(t) = : In addition, " # Z Z Z Z s s 2 s s E n n du = E n du + V dW n du s u u u u u 0 0 0 0 = x(s) + y(s) Thus, x(t) = 2 ( x(s) + y(s))ds and x (t) = 2 x(t) + 2 y(t): 33 Therefore, 2 2 t t 1 1 e + 4e 3 + 2 t x(t) = and the desired result is obtained. F Proof of Proposition 5 As dR = dS t t and dR = 2R dR + d[R] t t t 2(S S ) 1 t 0 = dS + d[S] ; t t 2 2 S S 0 0 we have 2(S S ) 2(S S )V t 0 t 0 t d[R; R ] = d[S] = dt: t t 3 3 S S 0 0 Note that Z Z Z Z t t t t p p s v E[S V ] = E S + n ds + V dW V + ( V )ds + V dW t t 0 s s 0 2 s s s s 0 0 0 0 Z Z t t = E S V + S ( V )ds + S V dW 0 0 0 2 s 0 s 0 0 Z Z Z Z Z Z t t t t t t + V n ds + t n ds n ds V ds + n ds V dW 0 s 2 s s 2 s s s 0 0 0 0 0 0 Z Z Z Z t t t t p p p s s s + V V dW + t V dW V ds V dW 0 s 2 s 2 s s s s s 0 0 0 0 Z Z t t p p s v + V dW V dW s s s s 0 0 Z Z Z Z t t t t = E S n ds V ds + n ds V dW 0 s 2 s s s 0 0 0 0 Z Z Z Z t t t t p p p s s v V ds V dW + V dW V dW : 2 s s s s s s s 0 0 0 0 For the last equality, the following assumption is used: E[V ] = V = ; E[n ] = 0: s 0 s The following can be derived: Z Z Z Z t t t s p p s s w(t) := E V ds V dW = E V V dW ds s s s u s u 0 0 0 0 Z Z Z Z Z t s s s s p p p s v s = E ( V )du V dW + V dW V dW ds 2 u u u u u u u 0 0 0 0 0 Z Z Z Z Z t s s t s = E V du V dW ds + E[V ]duds 2 u u u 0 0 0 0 0 = ( w(s) + s)ds = ( t 1 + e ): 34 Similarly Z Z Z Z t t t s p p v v q(t) := E n ds V dW = E n V dW ds s s s u s u 0 0 0 0 Z Z Z Z Z Z t s s t s s p p p v v s = E n du V dW ds + E V dW V dW ds 1 u u u u u u u 0 0 0 0 0 0 = ( q(s) + s)ds = ( t 1 + e ): Let Z Z Z Z Z Z t t t s t s z(t) := E n ds V ds = E n du V ds + E V du n ds s s u s u s 0 0 0 0 0 0 and Z Z Z Z s s s s E V n du = E V + ( V )du + V dW n du s u 0 2 u u u 0 0 0 0 Z Z Z Z s s s s = E V du n du + V dW n du 2 u u u u 0 0 0 0 = z(s) + ( s 1 + e ) 2 1 and Z Z Z Z s s s s E n V du = E n du + V dW V du s u u u u 0 0 0 0 = z(s) + ( s 1 + e ): Thus, s s z(t) = ( + )z(s) + (s 1 + e ) + ( s 1 + e ) ds 1 2 2 0 2 and 2 2 2 2 2 t 2 t ( + )t 1 2 1 2 + ( + )t + ( + )e + ( + )e e 1 2 2 1 1 2 1 2 1 2 1 2 1 2 2 1 z(t) = : 2 2 ( + ) 1 2 1 2 Therefore, t t 1 2 E[S V ] = S z(t) + ( t 1 + e ) ( t 1 + e ) + t t t 0 2 1 2 and hence t 2 2 1 2 2 2 t 2 t 1 2 E[S V ]du = S t + t + t t + 1 e t t + 1 e u u 0 1 2 2 2 2 0 2 2 2 2 2 2 ( )t + ( + )t 1 2 2 1 1 2 1 2 ( + ) 2 2 1 2 2 2 + + 1 2 1 2 1 2 2 t 1 t ( + )t 1 2 1 2 (e 1) (e 1) + (e 1) : 1 2 1 2 Finally, 2 2t E[[R; R ] ] = E[S V ]du : t u u 3 2 S S 0 0 0
Journal
Quantitative Finance – arXiv (Cornell University)
Published: Aug 14, 2019