Investigating the Statistical Properties of the Hurst Exponent Estimator of Rough Volatility Model Department of Finance Master’s Thesis in Finance Author: Advisor: Saeedeh Ostovari Dr. Alexander Herbertsson May 2021 Abstract The aim of this thesis is to provide a characterization of the statistical proper- ties of estimator of the Hurst parameter of the rough stochastic volatility model following fractional Brownian motion with Hurst index H. For this purpose, we perform a simulation experiment for fractional Brownian motion based on the cir- culant embedding method. Moreover, the study contributes to make a comparison between the Hurst estimator and the memory parameter estimator, d. The results indicate that the Hurst estimator is superior to considered memory es- timators, however, in the presence of microstructure noise, it is downward biased. Keywords: fractional Brownian motion, rough stochastic volatility models, cir- culant embedding method, fractionally integrated process, Realized volatility. 1 Acknowledgement Firstly, I would like to express my sincere gratitude to my advisors Alexander Herbertsson from Gothenburg University and Tommaso Proietti from Tor Vergata University for their continuous help and support. In addition, I am thankful to my teachers at the both Universities for their endeavors. I owe my deep gratitude to my siblings and my parents for their unconditional love, kindness, and support. Last but not least, my heartfelt gratitude has been reserved to my husband, Aboozar, who has been consistently supportive, without your love and patience this would not be possible to achieve. 2 Contents Abstract 6 Acknowledgment 6 1 Introduction 7 2 Literature Review 10 2.1 Rough Volatility . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.2 Fractionally integrated process . . . . . . . . . . . . . . . . . . . . 12 3 Methodology 13 3.1 introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 3.2 Fractional Brownian Motion . . . . . . . . . . . . . . . . . . . . . . 13 3.2.1 Self-Similarity . . . . . . . . . . . . . . . . . . . . . . . . . . 14 3.2.2 Dependent and Stationary Increments . . . . . . . . . . . . 14 3.2.3 Continuity . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.2.4 Fractional Gaussian Noise . . . . . . . . . . . . . . . . . . . 16 3.3 Rough Fractional Stochastic Volatility . . . . . . . . . . . . . . . . 17 3.3.1 Fractional Volatility . . . . . . . . . . . . . . . . . . . . . . 17 3.3.2 Scaling Property of Empirical RV . . . . . . . . . . . . . . . 18 3.4 Memory Behaviour . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 3.5 Simulation of fBm . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4 Estimation of Hurst Parameter 29 4.1 Monte Carlo Simulation . . . . . . . . . . . . . . . . . . . . . . . . 29 4.2 Estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 3 CONTENTS CONTENTS 4.2.1 Absolute Moments method . . . . . . . . . . . . . . . . . . . 33 4.2.2 Periodogram methods . . . . . . . . . . . . . . . . . . . . . 33 4.2.3 Statistical Property . . . . . . . . . . . . . . . . . . . . . . . 34 4.2.4 Estimation Analysis . . . . . . . . . . . . . . . . . . . . . . 35 5 Conclusion 40 Reference 41 Appendices 43 4 List of Figures 3.1 Paths of fBm for different values of H . . . . . . . . . . . . . . . . . 16 3.2 logm(q,∆) as a function of log ∆ . . . . . . . . . . . . . . . . . . . 20 3.3 Scaling of ζq with q, the slop is estimated H . . . . . . . . . . . . . 21 3.4 Histograms of log σt+∆− log σt for various Lags ∆; normal fit in black 21 3.5 The spectral density of fractional Gaussian noise process for differ- ent values of H. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 4.1 Mon(te ca(r∑lo of ρ̂(5)))as an estimator of ρ(5) . . . . . . . . . . . . . . 30 4.2 log Var m Hj=o bj is proportional to 2H . . . . . . . . . . . . . . 31 4.3 Scaling property of fractional Brownian motion with H = 0.1. The slope is close to 2H, 0.21. . . . . . . . . . . . . . . . . . . . . . . . 32 4.4 Scaling property of fractional Brownian motion contaminated with microstructure noise. . . . . . . . . . . . . . . . . . . . . . . . . . . 38 1 Paths of fGn for different values of H . . . . . . . . . . . . . . . . . 45 5 List of Tables 3.1 Estimated Hurst parameter H for log(RV ) series of 31 indices . . . 19 4.1 Slope of log-log plot of variance of partial sum versus sample size . 31 4.2 The statistical measures for the estimated Hurst parameter via ab- solute moments method . . . . . . . . . . . . . . . . . . . . . . . . 36 4.3 The statistical measures for the memory parameter d = H−0.5 via log-periodogram . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 4.4 The statistical measures for the memory parameter d = H−0.5 via GPH method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 4.5 The statistical measures for the memory parameter d = H−0.5 via GPH method, R = N . . . . . . . . . . . . . . . . . . . . . . . . . . 37 2 4.6 Estimated Hurst index H in presence of microstructure noise. . . . 39 4.7 Signal to Noise Ratio . . . . . . . . . . . . . . . . . . . . . . . . . . 39 6 Chapter 1 Introduction Volatility is a feature of financial markets. It is widely accepted to consider the financial asset return volatility as a stochastic process rather than a deterministic one. The process is not observable, so discovering the properties is not easy. A con- siderable empirical problem in pricing financial derivatives and risk management is how to extract and predict the volatility. Appropriate modelling of volatility is of foremost importance for financial applications. This can be applied to provide accurate interval predictions and to enhance density forecasts. Moreover, proper volatility prediction can be used for portfolio management. Fractional Brownian motion models have been applied to a wide range of ap- plications, for example, option pricing, hedging financial derivatives, and modeling volatility. Thus, it is important to provide a reliable estimator for its key param- eter, the fractal index, and make inference about the parameter. A seminal paper " Volatility is Rough", by Gatheral, Jaisson, and Rosenbaum (2018), provides empirical evidence that log-realized volatility behaves as fractional Brownian motion with Hurst parameter of order 0.1 that is called Rough fractional stochastic volatility, RFSV. The main focus of my thesis is to evaluates the statistical properties of the es- timator of the Hurst parameter of rough volatility proposed by Gatheral, Jaisson, and Rosenbaum (2018) by performing a simulation experiment. For this purpose, I applied the circulant embedding method, which is asymptotically exact Davies and Harte method in Davies and Harte (1987). In this method, the covariance 7 CHAPTER 1. INTRODUCTION matrix of Fractional Gaussian noise process is embedded in a circulant matrix in order to generate the fractional Gaussian noise, fGn, process. In fact, fractional Brownian motion is obtained by cumulating the fGn process. The fractional Brownian motion is a special case of Brownian semi-stationary pro- cess in which the Hurst index controls both roughness and long memory behaviour that cannot be disentangled Bennedsen, Lunde, and Pakkanen (2021). Proietti (2016) investigates modelling and predicting realized volatility series from a different point of view by proposing a new integrated moving average model in which the lag operator, B,1 is replaced by a fractional lag operator, 1 − (1 − B)d, called FLagIMA. The differencing degree known as a memory parameter is indicated by d. The parameter d is estimated by a Whittle-likelihood estimator which is consistent and asymptotically normal. In addition, this study aims to provide an in-deep analysis to make a com- parison between the Hurst estimator and that of memory parameter d. We apply periodogram methods to estimate the memory parameter, although there are other parametric and non-parametric methods. In this thesis, we find that the Hurst estimator applied in Gatheral, Jaisson, and Rosenbaum (2018) outperforms the considered memory estimators. The Hurst estimator is asymptotically unbiased when it is less than half, while it is downward biased in the presence of microstructure noise, the bias amplifies as the variance of microstructure noise increases. The main contribution of the study is to understand the statistical properties of the Hurst estimator with or without the presence of microstructure noise. As the volatility prediction is applied in, for example, option pricing and risk management, this finding may help to realize more precise volatility prediction in which volatility is driven by fractional Brownian motion. The rest of this study is structured as follows. Chapter 2 reviews the literature. It is divided in 2 subsections: Rough Volatility and fractionally integrated process. Chapter 3 provides methodology applied in this study divided in 5 parts: properties and definition of fractional Brownian motion, rough stochastic volatility model, memory behaviour, and simulation method; moreover, the result of testing the 1B is a backward operator such that (1−B)Yt = Yt − Yt−1 8 CHAPTER 1. INTRODUCTION scaling property for Realized Volatility series of 31 stock indices provided by the Oxford-Man Institute’s "realized library" are presented in this chapter. Chapter 4 gives analyses of estimators, and chapter 5 provides conclusions. 9 Chapter 2 Literature Review This chapter provides a review of some studies on modelling realized volatility. 2.1 Rough Volatility Daily data is commonly used to model volatility. The concept of realized volatil- ity is used when high frequency intraday data are applied to forecast and model volatility. Realized volatility, RV, is a daily non-parametric volatility measure ob- tained as a sum of squared intraday returns. To estimate the realized volatility as a proxy for variance at trading day t, assume that the day is split into m equally spaced intraday intervals, for example every 5 minutes, then the last observed price in each interval is considered as a closed price. The logarithmic price series {Pit}mi=0 is obtained by taking the logarithm of the closed prices. It is worth men- tioning that, as the prices are tick-by-tick transaction data, the series {P mit}i=0 can be formed by price interpolation. The realized volatility (RV) is calculated as the sum of intraday high frequency squared returns, ∑m (m) RV 2t = (Pit − Pi−1,t) . i=1 This estimator consistently estimates quadratic variation when the prices are ob- served without noise. To alleviate the impact of micro-structure noise, the high frequency data are sampled at low frequency, for example every 1-5 minutes. The 10 2.1. ROUGH VOLATILITY CHAPTER 2. LITERATURE REVIEW larger interval discovers more information of the true price volatility. The Oxford- Man Institute provides the daily realized variance as a proxy for daily volatility. Andersen and Bollerslev (1997) and Ding, Granger, and Engle (1993) provide evidence for long memory property of volatility process. The feature has been accepted as a stylized fact. Many researches have been conducted on modeling and forecasting volatility based on the assumption that long range dependence is a characteristic of volatility in the sense that autocorrelation function, ACF, is not integrable, see Andersen, Bollerslev, et al. (2003). Comte and Renault (1998) inspired by the empirical evidence of persistence in volatility discovered a continuous-time stochastic volatility model with the long memory property that follows a fractional Brownian motion process with Hurst parameter H > 0.5, which is called Fractional Stochastic Volatility, FSV. The notion of rough volatility has been introduced in a paper by Gatheral, Jais- son, and Rosenbaum (2018). The empirical studies have found that log realized variance time series of 21 indices from Oxford-Man institute’s Realized Library dataset are driven by a fractional Brownian motion, fBm, with H ≈ 0.1. They apply the fractional Brownian motion to forecast and model log realized volatility time series. They are inspired by the scaling property and stationarity of incre- ments of the fBm to estimate the Hurst index. Their finding provides a stark contrast to the stylized fact that volatility is a long memory process. Fukasawa, Takabatake, and Westphal (2019) investigate the roughness of realized volatility. They prove that the proposed quasi-likelihood estimator is consistent under high frequency asymptotics and apply it to time series of realized volatility. Akin to Gatheral, Jaisson, and Rosenbaum (2018), they obtain the same result, with a smaller rough parameter than that considered in Gatheral, Jaisson, and Rosenbaum (2018). Bennedsen, Lunde, and Pakkanen (2021) are inspired by Gatheral, Jaisson, and Rosenbaum (2018) and extend the analysis to a wider range of assets and draw the conclusion that volatility is rough and has heavy tails. They also introduce a new class of models so-called Brownian semi-stationary process for volatility of asset prices that can capture both roughness and strong persistence. Moreover, the models decouple roughness properties and long memory. The Hurst index H 11 2.2. FRACTIONALLY INTEGRATED PROCESS CHAPTER 2. LITERATURE REVIEW of fractional Ornstein-Uhlenbeck process can control both roughness and the rate of ACF decay when H ∈ (0, 0.5). 2.2 Fractionally integrated process The strong persistence in the auto-correlation function of Realized Volatility is a widespread stylized fact. Some studies have been conducted on modeling RV as a fractionally integrated process. A process {Yt} is a fractionally integrated process, also called ARFIMA(p,d,q)1, if ∆dYt is a ARMA(p,q)2 process, that is a fractionally integrated process of order d, I(d). ARFIMA(p,d,q) process is an extended version of ARIMA(p,d,q) process with a fractional value of d, instead of an integer value. This is an ordinary Autoregressive Moving Average process if d = 0. This kind of processes are strongly persistent and present long memory whenever the memory parameter d is positive and 0 < d < 1 , while they have 2 anti-persistence for −1 < d < 0, in Ruppert and Matteson (2015). 2 Baillie et al. (2019) investigate long memory property of RV series by applying fractionally integrated long memory models. They discover that the long memory processes can provide statistical interpretation of the persistence in RV. Proietti (2016) investigates modelling and predicting realized volatility series from a different point of view. He proposes a new integrated moving average model that the lag operator, B, is replaced by a fractional lag operator, 1−(1−B)d, called FLagIMA. The long memory parameter d is estimated by a Whittle-likelihood es- timator which is consistent and asymptotically normal. In the article, the realized volatility process is decomposed into two components, white noise, and a summa- tion of a fractional noise, through a generalized Beveridge-Nelson decomposition. At the end, he considers a collection of models to assess the forecasting perfor- mance of the FLagIMA using a Model Confidence Set (MCS). One of the contribution of this master study is to compare the memory parameter and the roughness parameter. 1It is defined in chapter 3. 2It is defined in the Appendix. 12 Chapter 3 Methodology 3.1 introduction In this chapter, we are going to establish the rough volatility in Gatheral, Jaisson, and Rosenbaum (2018) is driven by fractional Brownian motion, fBm. Therefore, the first section of this chapter defines the fBm and describes some properties ap- plied to estimate the rough index. The second part is dedicated to rough fractional stochastic volatility. The final part of the chapter presents one simulation method so-called circulant embedding method that is employed to assess the statistical properties of the roughness estimator. The contents in this chapter come from Dieker (2004), Shevchenko (2015), Mandelbrot and Van Ness (1968), Ruppert and Matteson (2015), Giraitis, Koul, and Surgailis (2012), and Brockwell and Davis (2016). 3.2 Fractional Brownian Motion Definition 3.1. A fractional Brownian motion process BH = {BH(t) : 0 ≤ t <∞} with H ∈ (0, 1) is a centered and self-similar Gaussian process character- ized by the following properties. 1)∀t ≥ 0, BHt has stationary increments. 2)∀t ≥ 0, BH0 = 0 and EBHt = 0, where symbol E denotes expectation. 3)∀t ≥ 0, EB2Ht = t2H . 13 3.2. FRACTIONAL BROWNIAN MOTION CHAPTER 3. METHODOLOGY 4)∀t > 0, has a Gaussian distribution. The covariance function is 1( E(BHBH) = t2H + s2Ht s − |t− s|2H). (3.1)2 The H parameter is called Hurst index or Hurst parameter. For a fixed value of H, the mean and covariance structure uniquely specify the distribution of the Gaussian process. According to the definition, for H = 1 the covariance is 2 E(BH H Ht Bs ) = t ∧ s, which that means Bt is regular Brownian motion process. Therefore, for H 6= 1 , the process is a generalization of Brownian motion so-called 2 fractional Brownian motion. Some properties are deduced from the definition (3.1). 3.2.1 Self-Similarity Definition 3.2. A homogeneous function g(x, y) of degree c > 0 is a real- valued function that for some constant c and ∀t ∈ R satisfies g(tx, ty) = tcg(x, y). (3.2) It is easy to check from the property (3.2) that the degree of homogeneity of the covariance function of fBm is 2H. Self-similarity. Now it can be clearly concluded from equation (3.1) that {BH(at) : 0 ≤ t <∞} and {aHBH(t) : 0 ≤ t <∞} for a > 0 have the same covariance and finite-dimensional distributions. That is called self-similarity of fBm with Hurst index H that means fBm is scale invariant. 3.2.2 Dependent and Stationary Increments The covariance function of fBm increments can be obtained from (3.1) and repre- sen[ted as follows ( E (BH−BH)(BH H 1t s t −Bs )] = |t1−s |2H2 +|t −s |2H−|t −t |2H2 1 2 1 −|s2−s |2H1 ) (3.3)1 1 2 2 2 Stationary Increments. The process Xt = BH Ht+s − Bs for t ≥ 0, based on (3.3) has the same covariance function as that of BHt . As for Gaussian process mean and 14 3.2. FRACTIONAL BROWNIAN MOTION CHAPTER 3. METHODOLOGY covariance determine the distribution, so the process has the same distribution as BHt . Dependent Increments. Let us look at the formula (3.3) and assume that s1 < t1 < s2 < t2 such that [s1, t1] ∩ [s2, t2] = ∅, no intersectio(n. Th)erefore,[ E (BHt −BHs1)(BHt −BH 1 s2)] < 0 for H ∈ 0,[ 1 2 ( 2)1 E (BHt −BHs )(BH1 t −BHs2)] > 0 for H ∈ , 11 2 2 Thus, the fBm is anti-persistent for H ∈ (0, 1), meaning that if the process is 2 decreasing, it would be more possible to increase in the future and vice versa, while in contrast, the fBm is persistent for H ∈ (1 , 1), refers to the long-range 2 dependence property of fBm, it means that the process tends to keep its trend in the future than to change it. The increments satisfy for any S > 0, t ∈ R, q > 0: E[|BH H q qHt+S −Bt | ] = κqS , where κq is the q-th moment of the absolute standard normal variable. The power of S is a linear function with the slope H, Therefore, the increments enjoy a considerable scaling property. 3.2.3 Continuity The continuity of fractional Brownian motion can be established based on the following formula, E[(BH −BH)2] = |t− s|2Ht s . (3.4) Theorem 3.1. Kolmogorov-Chentsov Continuity Theorem. Assume that for a stochastic process {Xt, t ≥ 0} there exist such K > 0, p > 0, β > 0 such that for all t ≥ 0, s ≥ 0 E[|Xt −Xs|p] ≤ K|t− s|1+β then the process X has a continuous modification, i.e a process {X̃t, t ≥ 0} such that X̃ ∈ C[0,∞) and for all t ≥ 0 Pr(Xt = X̃t) = 1. Moreover, for any γ ∈ (0, β )p 15 3.2. FRACTIONAL BROWNIAN MOTION CHAPTER 3. METHODOLOGY and and T > 0 the process X̃ is γ-Hölder continuous on [0, T ], i.e. |X̃t − X̃s| sup0≤s 1 . Referring to the article Comte and Renault 2 (1998), FSV is obtained by considering a truncated version of the general fBm. 17 3.3. ROUGH FRACTIONAL STOCHASTIC VOLATILITY CHAPTER 3. METHODOLOGY Assume that {St} is the price series of a given asset that is formulated as dSt = σtdωt St d log(σt) = νdω H t + α(µ− log(σt))dt E[dωHt , dωt] = ρdt, where σt is volatility, ν is the volatility of volatility, µ and α are mean reversion parameters. ωt is ordinary Brownian motion and ωHt is a truncated version of fBm with H ∈ (1 , 1) that is defined by 2 ∫ t (t− s)H− 12 ωH(t) = Γ(H + 1 dω(s), 0 )2 where Γ is the Gamma function. The process is called type II fractional Brownian motion. Marinucci and Robinson (1999) assessed that the incremental process of the type II fBm is non-stationary, while increments of the type I fBm are stationary, the process is defined∫ast − H− 1 ∫H (t s) 02 (t− 1 1s)H− 2 − (−s)H− 2ω1 (t) = dω(s) + dω(s). 0 Γ(H + 1) −∞ Γ(H + 1) 2 2 3.3.2 Scaling Property of Empirical RV Gatheral, Jaisson, and Rosenbaum (2018) take advantage of the scaling property of fBm to estimate the Hurst index. For this purpose, they define a quantity m(q,∆) as: 1 ∑n m(q,∆) = | log(σk∆)− log(σ q(k−1)∆)| (3.7) n k=1 that {σt} is a daily realized volatility series on a time grid with mesh ∆ on interval [0, T ], such that k ∈ {1, 2, ..., b T c}, n = b T c and q ≥ 0. Under the assumption ∆ ∆ that for some bq > 0 and sq ≥ 0, the following holds N qsqm(q,∆)→ bq , ∆→ 0 (3.8) the smoothness parameter sq controls regularity of the volatility process. To estimate sq, it is needed to calculate m(q,∆) for different value of ∆ and q, 18 3.3. ROUGH FRACTIONAL STOCHASTIC VOLATILITY CHAPTER 3. METHODOLOGY then run a linear regression of log m(q,∆) versus log ∆. As the volatility is not observable, we use 5-minute realized volatility series as a proxy for the spot volatility. The RV series of 31 indices from January 3, 2000 to November 2020 are taken from Oxford-Man Institute of Quantitative Finance Realised Library 1. Figure 3.2 exhibits the scaling property of RV series of FTSE and S&P indices as the plots of log m(q,∆) against log ∆ for different values of q which can be described as E[| log(σ∆+t)− log(σt)|q] = b ∆ζqq , ∀t > 0. where ζq = qsq is the slope of the regression and a linear function of q. The estimation of H is the coefficient parameter in regression of ζq against q. Table 3.1 exhibits the estimated values of Hurst exponent for realized volatility series of 31 indices. Table 3.1: Estimated Hurst parameter H for log(RV ) series of 31 indices index Ĥ index Ĥ index Ĥ index Ĥ FTSE 0.11 N225 0.11 DAX 0.12 SPX 0.13 RUT 0.09 AEX 0.13 AORD 0.093 BFX 0.13 BSESN 0.10 BVSP 0.13 DJI 0.13 FCHI 0.12 HSI 0.09 IBEX 0.11 IXIC 0.13 KS11 0.12 KSE 0.08 MXX 0.08 NSEI 0.09 SSMI 0.15 STI 0.09 STOX 0.06 BVLG 0.14 FTMIB 0.12 GSPTSE 0.11 SSEC 0.12 OMXC20 0.09 OMXHPI 0.10 OMXSPI 0.11 OSEAX 0.11 SMSI 0.10 - - The figure 3.3 represents the linear relationship that is called mono-fractal scaling property. Moreover, Gatheral, Jaisson, and Rosenbaum (2018) provides empirical evidence that the increments of logarithmic realized volatility have nor- mal distribution by fitting a normal distribution to the empirical data, figure 3.4 visualises it. Therefore, based on the two properties of log realized volatility the 1https://realized.oxford-man.ox.ac.uk/data/download 19 3.3. ROUGH FRACTIONAL STOCHASTIC VOLATILITY CHAPTER 3. METHODOLOGY incremental process may be constructed as the following log σ H Ht+∆ − log σt = ν(Bt+∆ −Bt ) that BHt is fBm with Husrt index H. The above equation can be written as σt = σ exp(νB H t ) where σ is positive. Therefore, log σt has the same distribution as fractional Brow- nian motion. Indeed, Gatheral, Jaisson, and Rosenbaum (2018) apply absolute moment methods to estimate the Hurst parameter via the incremental process by employing the stationarity and scaling property of increments. Figure 3.2: logm(q,∆) as a function of log ∆ 20 3.4. MEMORY BEHAVIOUR CHAPTER 3. METHODOLOGY Figure 3.3: Scaling of ζq with q, the slop is estimated H Figure 3.4: Histograms of log σt+∆ − log σt for various Lags ∆; normal fit in black 3.4 Memory Behaviour The memory property of a stationary process can be described in terms of depen- dence structure or correlation across time series. The long-range dependence or 21 3.4. MEMORY BEHAVIOUR CHAPTER 3. METHODOLOGY long memory refers to strong correlation in a time series. The dependence can be measured via spectral density in frequency domain or the autocovariance function in time domain. In other words, time-domain analysis of long and short memory relies on the asymptotic behaviour of autocovariance function γ(k), as k → ∞, while frequency domain analysis relies on the periodogram or sample spectral den- sity, which is discrete Fourier transform of the autocovariance sequence. Definition 3.4. (Spectral density) Suppose that {Yt} is a stationary process with autocovariance γ(.). The spectral density of {Yt} is the Fourier transform of the autocovariane sequence defined as 1 ∑∞ f(ω) [= γ(j) exp (−iωj)2π j=−∞∑ ]∞1 = γ(0) + 2 γ(j)cos(ωj) 2π j=1 √ where eiω = cos(ω)+i sin(ω) and i = −1, ω ∈ (−π, π]. This is a non-negative and periodic function with period 2π, also it is an even function in Brockwell and Davis (2016). Definition 3.5. A stationary process {Yt} with covariance function γ(k) has long memory or long-range depen∑dence if |γ(k)| =∞; (3.9) k∈Z short memory if ∑ ∑ |γ(k)| <∞ and γ(k) > 0; (3.10) k∈Z k∈Z anti-persistence, or neg∑ative memory if ∑ |γ(k)| <∞ and γ(k) = 0. (3.11) k∈Z k∈Z The summability condition (3.10) implies that the spectral density f is bounded, while under the anti-persistence condition (3.11) f(0) = 0. For a long memory process, the spectral density is not bounded at the zero frequency. 22 3.4. MEMORY BEHAVIOUR CHAPTER 3. METHODOLOGY The memory property also can be characterised by the asymptotic behaviour of autocovariance function that hyperbolically decays to zero as follows, γ(k) = cγ|k|−1+2d 1 , ∀k > 1, ∃ 0 < d < , cγ > 0. (3.12) 2 Threfore, if the power 1−2d in (3.12) lies in the interval (0, 1), the pr∑ocess has long memory, means that 0 < d < 1 . Assume that −1 < d < 0 and either k∈Z γ(k) = 02 2 or γ(k) = −|k|−1+2dcγ, then γ(k) holds the negative memory condition 3.11. Definition 3.6. (Memory parameter). The parameter d in (3.12) is called memory parameter of a stationary process that is positive for a long memory, and negative for a negative memory process. The parameter is zero if the stationary process has short memory. The definition below provides memory concepts in the frequency domain. Definition 3.7. Suppose a stationary process has a spectral density that is bounded on [ζ, π] for any ζ > 0, and satisfies f(ω) ∼ c |ω|−2d −1 1f as ω → 0 for < d < , cf > 0. (3.13) 2 2 The memory characteristic of a stationary process depends on whether d = 0, d ∈ (−1 , 0), or d ∈ (0, 1) has short memory, negative memory, or long memory 2 2 respectively. As it was mentioned in subsection 3.2.2, fractional Brownian motion with Hurst parameter H ∈ (0.5, 1) has long range dependence. The process has stationary increments that are anti-persistent for H ∈ (0, 0.5). According to the equations (3.6) and (3.12), Hurst parameter H and memory parameter d are related such that H = d+ 1 . 2 The figure 3.5 below can visualize the memory property of incremental process of fBm in the sense of spectral density behaviour around zero frequency, the right side plot exhibits the long-range dependence, and the left side one represents the anti-persistence of the increments of fBm process, or fGn. 23 3.4. MEMORY BEHAVIOUR CHAPTER 3. METHODOLOGY (a) H < 0.5 (b) H > 0.5 Figure 3.5: The spectral density of fractional Gaussian noise process for different values of H. Definition 3.8. The class of Fractionally integrated process, ARFIMA(p,d,q), is an extended version of ARIMA(p,d,q) process by the difference that the d is fractional. Consider the stochastic process φ(B)(1−B)d(Yt − µ) = θ(B)t, where d is a real parameter and φ(z) = 0 and θ(z) = 0 for | z |> 1. {t} is a white noise with mean zero and variance σ2, and B is the lag operator. • The process is nonstationary if d ≥ 0.5. • It is stationary and invertible for −1 < d < 1 , in which case we say that Y 2 2 t is a fractionally integrated ARMA process, ARIMA(p,d,q). • The process is said to be persistence if 0 < d < 0.5 and anti-persistence if −0.5 < d < 0. The ARFIMA(0,d,0) is called Fractional noise process represented as ∑∞ Γ(j + d) Yt = (1−B)−d 2t = t−j, t ∼ WN(0, σ ). (3.14) Γ(d)Γ(j + 1) j=0 24 3.5. SIMULATION OF FBM CHAPTER 3. METHODOLOGY Γ(x) is the Gamma function. The autocorrelation is 2 Γ(h+ d)Γ(1− d)ρ(h) = σ (3.15) Γ(d)Γ(1 + h− d) The spectral density is obtained by the squared gain of the filter (1−B)−d times the spectral density of the white noise σ2t, f = ,2π f(ω) = f | 1− e−iω |−2d = f[2(1− cosω)]−d (3.16) | ω= f 2 sin( ) |−2d 2 The expression | 2sin(ω ) |−2d converges to |ω|−2d when ω → 0, then we get the 2 equation 3.13. 3.5 Simulation of fBm In the previous sections, we have seen properties of fractional Brownian motion that is a self-similar and Gaussian process applied to model Realized Volatility series. In this part, we aim to elaborate how the process is simulated in order to assess the statistical property of the Hurst estimator. The process only can be simulated in discrete time, so the discrete value of fractional Brownian motion is denoted by {BH Nn }n=1 that N is sample size. A realization of fBm on interval [0, T ] for some T > 0 can be obtained via applying the self-similarity property ,i.e. by multiplying {BH}N T Hn n=1 by ( ) .N Davies and Harte (1987) proposed a method that was generalized by Andrew and Chan (1994) the so-called circulant embedding method. The idea behind the method is to exploit the Gaussian property of fGn. Assume that X is a Gaussian vector with specific mean µ and covariance Γ it can be obtained by a standard Gaussian vector Y as X = µ + SY where the matrix Γ can be decomposed such as Γ = SSᵀ. Therefore, the basic idea is to find the matrix S, the square root of the covariance matrix of fGn. 25 3.5. SIMULATION OF FBM CHAPTER 3. METHODOLOGY It is needed to assume that sample size N is of power two, N = 2g for some g ∈ N. Furthermore, to calculate the square root matrix S, the idea is to embed the covariance matrix Γ in a circulant matrix C of size 2g+1. According to the definition of fractional Gaussian noise (3.3) the covariance matrix is a Teoplitz matrix that has[the following form1 γ(K) = |K − 1|2H − 2|K|2H + |K + 1|2H ] (3.17) 2    γ(0) γ(1) · · · γ(N − 1)  γ(1) γ(2) · · · γ(N − 2) Γ = .. . . .. . . . ...  γ(N − 1) γ(N − 2) · · · γ(0) The circulant covarianc{e matrix is defined as, γ(i) if i = 1, 2, ..., N − 1 ci = γ(M − i) if i = N,N + 2, ...,M − 1. such that c0 = γ(0) = 1, cN = 0 and M = 2(N − 1).  c0 c1 c2 · · · cM−2 cM−1 · · ·    cM−1 c0 c1 cM−3 cM−2  c C M−2 cM−1 c0 · · · cM−4 cM−3 =  .. .. .. . . . .. ..  . . . . . c2 c3 c4 · · · c0 c1  c1 c2 c3 · · · cM−1 c0 Note that the matrix is symmetric and positive definite, generally circulant matrix is not positive definite, row j can be formed by shifting the first row j − 1 places, then the eliminated elements are added on the left. The algorithm exploits the following theorem, whose proof is available in Dietrich and Newsam (1997). Theorem 3.2. The circulant matrix C can be decomposed as C = Y ΛY ∗, where matrix Y is unitary with conjugate transpose matrix Y ∗. Λ is the diagonal matrix 26 3.5. SIMULATION OF FBM CHAPTER 3. METHODOLOGY constituted by eigenvalues of matrix C that are positive and real and defined as 2∑N−1 jk cj exp(2πi ) for k = 0, 1, 2, ..., 2N − 1. (3.18) 2N j=0 where cj is the j + 1 element of the first row of the circulant matrix. The unitary matrix Y is defined as 1 − jkyjk = exp( 2πi ) for j, k = 0, 1, 2, ..., 2N − 1. 2π 2N Since C = Y ΛY ∗ and Y ∗ is unitary, the matrix S = Y Λ1/2Y ∗ can satisfy SS∗ = C. Consequently, S satisfies the desired property. To simulate the fBm process, we need to simulate Y Λ1/2Y ∗W , where W is a standard normal vector. It is done in the following steps: 1- Computation of eigenvalues with equation (3.18) via the Fast Fourier Trans- form, FFT. 2- Calculate V = Y ∗W . This is done through the two following steps: • Generate two standard normal variables V0,VN ; • First generate independent standard normal random variables (1), (2)Wj Wj , then 1 (1) (2) Vj = √ (Wj + iWj ) 2 1 (1) (2) V2N−j = √ (Wj − iWj ) 2 3- In this step Z = Y Λ1/2V is computed as 2 1 ∑N−1√ jk Zk = √ λjVj exp(−2πi ) for k = 0, 1, 2, ..., 2N − 1. 2N 2N j=0 27 3.5. SIMULATION OF FBM CHAPTER 3. METHODOLOGY The sequence (Z)k=2N−1K k√=0 is the Fourier transform of  λ √ k V0 if k = 0;2N√ λ k (1) (2)  (Wk + iWK ) if k = 1, .., N − 1;w 4Nk =  √ λ k V2N N if k = N ;λk (1) (2)(W 4N 2N−k + iW2N−K) if k = N + 1, .., 2N − 1; The best way is FFT to speed up the computations. This method is advantageous in sense of speed, for N sample points the number of computation is N log(N). Finally, we have a sample fractional Gaussian noise by considering the first N elements of Z. The fractional Brownian motion is obtained by taking cumulative sum of fGn. 28 Chapter 4 Estimation of Hurst Parameter The previous chapter briefly introduced RFSV, fBm process, and also one of the simulation methods of fractional Brownian motion. Moreover, the results of per- forming the estimation method in Gatheral, Jaisson, and Rosenbaum (2018) on logarithmic realized volatility series of 31 indices illustrated that the RV is rough with H ≈ 0.1. This chapter presents the result of the study and discuss the findings. The sta- tistical properties of the Hurst estimator was assessed by performing simulation of fBm. I also hope that the simulated samples satisfy some properties of the exact sample. The first part of the chapter investigates the simulation results and the second part is dedicated to estimation methods and their properties. 4.1 Monte Carlo Simulation In section 3.5, we elaborated the method we applied to simulate fractional Brown- ian motion. we want to verify that whether the simulations of fBm can satisfy some properties of fBm process. For this purpose, we generate 1000 sample paths of fBm of length N = 212 for different values of Hurst parameter H ∈ {0.1, 0.05, 0.6, 0.8}. The first test aims to check whether the empirical ACF, ρ̂(r), approximate the theoretical ACF, ρ(r), at lag r. The figure 4.1 visualizes it at lag r = 5, means that the model is capable of approximating the true ACF structure of fractional 29 4.1. MONTE CARLO SIMULATION CHAPTER 4. ESTIMATION OF HURST PARAMETER Gaussian noise. It is worth mentioning that as fBm is cumulative sum of fGn, therefore the model can explain the true ACF structure of fBm. Figure 4.1: Monte carlo of ρ̂(5) as an estimator of ρ(5) The second test shows that the logarithmic variance of partial sum of the simulated sample paths of fGn is a linear function of logarithm of the sample size with slope equal to 2H. According to fGn definition 3.3, we have ∑m H H [ ( )bj∑j=0 ] = Bm m [ ( )] log Var bHj = log Var B H m (4.1) j=o = log(m2H) = 2Hlog(m) The figure 4.2 exhibits the linear relationship, the slopes is expected to be close 2H, Table 4.1 represents the slopes for different values of N ∈ {210, 212, 213}, it can be seen that as N increases the estimated slopes get closer to the true value 30 4.1. MONTE CARLO SIMULATION CHAPTER 4. ESTIMATION OF HURST PARAMETER of 2H. Table 4.1: Slope of log-log plot of variance of partial sum versus sample size N H = 0.1 H = 0.05 H = 0.6 H = 0.8 210 0.1897 0.0986 1.2182 1.6128 212 0.1909 0.0883 1.2013 1.6212 213 0.2150 0.0993 1.2133 1.6109 ( (∑ )) Figure 4.2: log Var m bHj=o j is proportional to 2H 31 4.2. ESTIMATION CHAPTER 4. ESTIMATION OF HURST PARAMETER The figure 4.3 depicts the scaling property when N = 213 and H = 0.1. There- fore, the three tests provide evidence that the simulation method works well. Figure 4.3: Scaling property of fractional Brownian motion with H = 0.1. The slope is close to 2H, 0.21. 4.2 Estimation The long range dependence can be modeled by fractionally integrated models, fBm, and fGn models. For fractional Brownian motion or fractional Gaussian noise, the Hurst exponent measures the dependence level, when 0.5 < H < 1 shows persis- tence. In the fractionally integrated process, ARFIMA(p,d,q), the persistency is measured via the differencing operator, d ∈ (−0.5, 0.5). As mentioned in Benned- sen, Lunde, and Pakkanen (2021), the Hurst index controls both the long-term and roughness behaviour of a process through H = 1 + d. The relationship between 2 the parameters can be achieved by putting equal the exponents in equations 3.6 and 3.10 so we have 2H − 2 = α = 2d− 1 → 1H = + d. (4.2) 2 This study applies the Absolute Moments method in Gatheral, Jaisson, and Rosen- baum (2018) and the Periodogram methods to estimate Hurst parameter H, and memory parameter d, then confirms the relationship, H = 1 + d by comparing the 2 32 4.2. ESTIMATION CHAPTER 4. ESTIMATION OF HURST PARAMETER estimation results. The performance of estimators are assessed by their statistical properties. 4.2.1 Absolute Moments method The absolute moments method applies the principle that the increments of fBm process, {BH − BHt s }s 0.5 or d > 0. Part 1: The bias of Hurst estimator is considerably smaller than bias of other estimators. The estimator is asymptotically unbiased as sample size increases. Moreover, its MSE is lower than GPH estimator and log-periodogram estimator. Therefore, it can be concluded that the absolute moments method outperforms the other es- timators when H < 0.5. The first two rows of table 4.4 shows that the GPH estimator may not be a good estimator for Hurst parameter in the region. Part 2: In case when H > 0.5, the bias of log-periodogram estimator is smallest, while the MSE is higher than others. As the three estimators are bias, the more effi- cient estimator is the one with lower mean square error. Both absolute moments estimator and GPH estimator has MSE less than 10−3 for N = 213. The GPH estimator depends on R, table 4.5 presents the statistical result of GPH estimator when R = N/2. The level of bias and MSE change as R changes. Comparing the two tables 4.5 and 4.2 shows that the MSE of the absolute moment estimator is smaller than the GPH estimator, therefore, the absolute moments estimator has the best performance. Table 4.5: The statistical measures for the memory parameter d = H − 0.5 via GPH method, R = N . 2 Mean Bias MSE N 210 212 213 210 212 213 210 212 213 d = −0.4 -0.5967 -0.5939 -0.5951 -0.1967 -0.1939 -0.1951 0.0398 0.0379 0.0382 d = −0.45 -0.7182 -0.7189 -0.7192 -0.2682 -0.2689 -0.2692 0.0732 0.0726 0.0726 d = 0.1 0.1249 0.1233 0.1226 0.0249 0.0233 0.0226 0.0017 0.0008 0.0006 d = 0.3 0.3623 0.3608 0.3600 0.0623 0.0608 0.0600 0.0050 0.0039 0.0037 It would be interesting to look at fractional Brownian motion with micro- structure noise and assess whether the new process has the scaling property and how biased is the Hurst estimator. For the purpose, it is supposed that the mi- crostructure noise (MN) has Gaussian distribution with mean 0 and variance σ2. 37 4.2. ESTIMATION CHAPTER 4. ESTIMATION OF HURST PARAMETER The MN is simulated independently from fBm, then added to fBm. The total variance of the new signal is fBm variance plus microstructure noise variance. Different variances are used to examine the Hurst estimator’s properties, such as how high the bias becomes as the noise variance grows. In the presence of mi- crostructure noise, the table (4.6) shows estimated Hurst parameter values. When σ = 0, the Hurst estimator has the lowest bias and MSE since the fractional Brow- nian motion is not contaminated by MN. It is obvious that the Hurst estimator is downward biased. As the variance of the signal or the variance of the MN increases, the biased and MSE increase. For contaminated fBm with the Hurst index H > 0.5, the result is the same. The signal-to-noise ratio, or the ratio of signal variance to noise variance, is one indicator of the amount of noise contamination in a signal. The signal becomes more contaminated with microstructure noise as the ratio decreases; table (4.7) represents the average ratio over 1000 signal replications. One of the properties of increments of fBm, as seen in the previous chapter, is that they have a scaling property. The increments of the new signal also have a scaling property, shown in the figure (4.4). Figure 4.4: Scaling property of fractional Brownian motion contaminated with microstructure noise. 38 4.2. ESTIMATION CHAPTER 4. ESTIMATION OF HURST PARAMETER Table 4.6: Estimated Hurst index H in presence of microstructure noise. Std Bias MSE σ1 = 0 -0.0023 0.0007 σ2 = 0.0005 -0.0483 0.0026 σ3 = 0.0006 -0.0649 0.0044 σ4 = 0.0008 -0.0757 0.0058 σ5 = 0.001 -0.0820 0.0068 σ6 = 0.0041 -0.0986 0.0097 Table 4.7: Signal to Noise Ratio Standard deviation σ2 σ3 σ4 σ5 σ6 Signal-to-Noise Ratio 2.7740 1.8796 1.5241 1.3494 1.0221 39 Chapter 5 Conclusion The study has aimed to investigate the statistical property of the Hurst estimator, which is used in Gatheral, Jaisson, and Rosenbaum (2018). Without microstruc- ture noise, the Hurst estimator outperforms the log-periodogram and GPH estima- tors, as we saw in the previous chapter. Furthermore, for H < 0.5, the estimator seems to be asymptotically unbiased. However, there is some bias in the presence of microstructure noise. The bias is downward, and the bias increases as the noise variance increases. As a conse- quence, as the noise variance grows, the Hurst estimator’s MSE increases. Gatheral, Jaisson, and Rosenbaum (2018) pointed out that the Hurst estimator is downward biased due to noise which is consistent with our findings. Understanding the characteristics of volatility is foremost important, since volatil- ity is an unobserved feature of financial markets and the volatility prediction is applied in portfolio management, risk management and pricing financial deriva- tives. 40 Bibliography Andersen, T.G. and T. Bollerslev (1997). “Intraday periodicity and volatility per- sistence in financial markets”. In: Journal of Empirical Finance 4.2, pp. 115– 158. Andersen, T.G., T. Bollerslev, et al. (2003). “Modeling and Forecasting Realized Volatility”. In: Econometrica 71.2, pp. 579–625. Andrew, T. A. Wood and G. Chan (1994). “Simulation of Stationary Gaussian Processes in [0,1]d”. In: Journal of Computational and Graphical Statistics 3.4, pp. 409–432. Baillie, R. et al. (2019). “Long Memory, Realized Volatility and Heterogeneous Autoregressive Models”. In: Journal of Time Series Analysis 40.4, pp. 609– 628. Bennedsen, M., A. Lunde, and M.S. Pakkanen (Jan. 2021). “Decoupling the Short- and Long-Term Behavior of Stochastic Volatility”. In: Journal of Financial Econometrics. Beran, J. (1994). Statistics for long-memory processes. Brockwell, p.J and R. Davis (2016). Introduction to Time Series and Forecasting. Springer texts in statistics. Springer. Comte, F. and E. Renault (1998). “Long memory in continuous-time stochastic volatility models”. In: Mathematical Finance 8.4, pp. 291–323. Corsi, F. (Feb. 2009). “A Simple Approximate Long-Memory Model of Realized Volatility”. In: Journal of Financial Econometrics 7.2, pp. 174–196. Davies, R. B. and D. S. Harte (Mar. 1987). “Tests for Hurst effect”. In: Biometrika 74.1, pp. 95–101. 41 BIBLIOGRAPHY BIBLIOGRAPHY Dieker, T. (2004). Simulation of fractional Brownian motion. MSc Thesis, Univer- sity of Twente, Amsterdam, The Netherlands. Dietrich, C. R. and G. N. Newsam (1997). “Fast and Exact Simulation of Stationary Gaussian Processes through Circulant Embedding of the Covariance Matrix”. In: SIAM Journal on Scientific Computing 18.4, pp. 1088–1107. Ding, Z., C. Granger, and R.F. Engle (1993). “A long memory property of stock market returns and a new model”. In: Journal of Empirical Finance 1.1, pp. 83– 106. Fukasawa, M., T. Takabatake, and R. Westphal (2019). “Is Volatility Rough”. In: arXiv: Statistics Theory. Gatheral, J., T. Jaisson, and M. Rosenbaum (2018). “Volatility is rough”. In: Quan- titative Finance 18.6, pp. 933–949. Geweke, J. and S. Porter-Hudak (1983). “The estimation and application of long memory time series models”. In: Journal of Time Series Analysis 4.4, pp. 221– 238. Giraitis, L., H. Koul, and D. Surgailis (2012). Large sample inference for long memory processes. Imperial College Press. isbn: 1848162782 9781848162785. Mandelbrot, B. and J. Van Ness (1968). “Fractional Brownian Motions, Fractional Noises and Applications”. In: SIAM Review 10.4, pp. 422–437. Marinucci, D. and P.M. Robinson (1999). “Alternative forms of fractional Brownian motion”. In: Journal of Statistical Planning and Inference 80.1, pp. 111–122. Proietti, T. (Apr. 2016). “Component-wise Representations of Long-memory Mod- els and Volatility Prediction”. In: Journal of Financial Econometrics 14.4, pp. 668–692. Ruppert, D. and D. Matteson (2015). Statistics and data analysis for financial engineering : with R examples. Springer texts in statistics. Springer. Shevchenko, G. (2015). “Fractional Brownian motion in a nutshell”. In: Interna- tional Journal of Modern Physics: Conference Series 36. 42 Appendices 43 Definition .1. (Weak Stationarity) A stochastic process {Yt}t≥0 is stationary (covariance stationary or weak stationary) if for ∀t • E(Yt) = µ <∞ • V(Yt) = σ2y <∞ • Cov(Yt, Yt+s) = γ(s) is not dependent on t. depends on t and the the cov(Xs, Xt) only depends on |t− s| for all t, s ∈ Z. Definition .2. (Strong Stationarity) A stochastic process {Yt}t≥0 is strongly stationary if its joint probability distribution does not change with shifting time, therefore its variance and mean stay unchanged over time. According to the definition .1 and .2, weak stationarity and Strong stationarity are equivalent given that the process is a Gaussian process as the distribution is determined only by its mean and variance. Definition .3. (ARMA) {Yt} is an autoregressive moving average process of degree p and, that is called ARMA(p,q) if it is stationary and if for every t, Yt − φ1Yt−1 − ...− φpYt−p = t + θ1t−1 + ...+ θqt−q, where {t} ∼ WN(0, σ2) and the polynomials (1 − φ1Z − ...φpZp) and (1 + θ1 + ...+ θqZ q) have no common factors. This can be written in a concise form as φ(B)Yt = θ(B)t, where φ(.) and θ(.) are the pth and qth-degree polynomials. The {Yt} is an autoregressive process of order p, AR(p), if θ(Z) ≡ 1, and a moving average process of order q, MA(q), if φ(z) ≡ 1, Brockwell and Davis (2016). Definition .4. (Teoplitz Matrix) A Teoplitz Matrix is a matrix in which each descending diagonal from left to right are constant. Any n × n matrix A matrix of the form bellow is called Teoplitz matrix if Aij = Ai+1,j+1 = ai−j 44 Aij is the i, j element of matrix A.  a0 a−1 a−2 · · · · · · a−(n−1) . . a1 a0 a−1 . a−(n−2)    .a a . . . . . . . ..  2 1 . . A =   ... . . . . . . . . .  a−1 a−2 ... . . . a1 a0 a−1  an−1 · · · · · · a2 a1 a0 Since the covariance function of a stationary process is symmetric, γ(−k) = γ(k), the Teoplitz matrix is symmetric. Figure 1: Paths of fGn for different values of H 45