{'fraction_non_alphanumeric': 0.04946477006853428, 'fraction_numerical': 0.02205105907199205, 'mean_word_length': 4.013806706114399, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 23, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Semiparametric inference on the fractal index of Gaussian and conditionally Gaussian time series data

We study a well-known estimator of the fractal index of a stochastic process. Our framework is very general and encompasses many models of interest; we show how to extend the theory of the estimator to a large class of non-Gaussian processes. Particular focus is on clarity and ease of implementation of the estimator and the associated asymptotic results, making it easy for practitioners to apply the methods. We additionally show how measurement noise in the observations will bias the estimator, potentially resulting in the practitioner erroneously finding evidence of fractal characteristics in a time series. We propose a new estimator which is robust to such noise and construct a formal hypothesis test for the presence of noise in the observations. Finally, the methods are illustrated on two empirical data sets; one of turbulent velocity flows and one of financial prices.


Introduction
Fractal-like models are used in a wide array of applications such as in the characterization of surface smoothness/roughness (Constantine and Hall, 1994), in the study of turbulence (Corcuera et al., 2013), and many others (e.g., Burrough, 1981;Mandelbrot, 1982;Falconer, 1990). Most recently, these models have attracted attention in mathematical finance, as models of stochastic volatility (e.g., Gatheral et al., 2018;Bayer et al., 2016;Bennedsen et al., 2017a,b;Jacquier et al., 2017). In such applications, it is imperative to be able to estimate and conduct inference on the key parameter in these models, the fractal index. Many estimators of this parameter exist (see Gneiting et al., 2012, for a survey); however, the underlying assumptions behind the various estimators, as well as their asymptotic properties, are often different and rarely stated in a clear and concise manner. These facts can make analysis difficult for the practitioner, as well as for the researcher.
This paper aims at making empirical analysis in applications, such as the ones mentioned above, easier. We clearly lay out a large and coherent framework -including the valid underlying assumptions -for analysing time series data which are potentially fractal-like. We focus on a specific estimator, which is arguably the most widely used in practice and which in our experience is the most accurate. Further, the estimator is easy to implement -it relies on a simple OLS regression -and its asymptotic properties are easy to apply. Our hope is that this will provide a transparent guide to analyzing fractal data using sound statistical methods.
The main contribution of the paper is to lay out the theory of the estimator and provide the theoretical underpinnings of it, stating the results in a manner so that application of the results becomes straight forward. For this, we rely heavily on earlier theoretical work on the increments of fractal processes, most notably  and Barndorff-Nielsen et al. (2011). We further investigate the estimator numerically to gauge it's properties when applied to data, leading to a number of practical recommendations for implementation. Most importantly, we advocate a different choice of bandwidth parameter for the estimator, than what is generally accepted practice in the literature, cf. Section 3.1.
In their survey of the asymptotic theory of various estimators of the fractal index, Gneiting et al. (2012) section 3.1., report that "a general non-Gaussian theory remains lacking". The second contribution of this paper is to extend the estimation theory beyond the Gaussian paradigm. We accomplish this by volatility modulation which turns out to be a convenient way of extending the theory to a large class of non-Gaussian processes. As will be seen, this results in conditionally Gaussian processes for which the fractal theory continues to hold. Again we clearly lay out the relevant assumptions and focus on the interpretation of the results and implementation of the methods.
The final contribution of the paper is an in-depth study of the case where the data are contaminated by noise, such as measurement noise. We prove that noise will bias estimates of the fractal index downwards, thereby making noise-contaminated data look more rough than the underlying process actually is. We go on to propose a novel way to construct an estimator which is robust to noise in the observations. The new estimator also relies on an OLS regression and is just as easy to implement as the standard (non-robust) estimator studied in the first part of the paper. We present the asymptotic theory concerning the robust estimator and propose a hypothesis test, which can be used to formally test for the presence of noise in the observations. The rest of the paper is structured as follows. Section 2 presents the mathematical setup and assumptions and gives some examples of the kind of processes we have in mind. The section then goes on to consider some extensions to the basic setup, most notably the extension to non-Gaussian processes. Section 3 presents the semiparametric estimator of the fractal index and it's asymptotic properties. Then, in Section 3.2, we consider the case where the observations have been contaminated by noise and present asymptotic theory for a new estimator in this case; Section 3.2.1 presents a formal test for the presence of noise. Section 4 contains small simulation studies, illustrating the finite sample properties of the asymptotic results presented in the paper. Finally, Section 5 contains two illustrations of the methods: the first using measurements of the longitudinal component of a turbulent velocity field, and the second using a time series of financial prices. Section 6 concludes and gives some directions for future study. Proofs of technical results and some mathematical derivations are given in an appendix.

Setup
Let (Ω, F, P) be a probability space satisfying the usual assumptions and supporting X, a onedimensional, zero-mean, stochastic process with stationary increments. Define the p'th order variogram of X: As we intend to make use of the theory developed in , 2011 we adopt the assumptions of those papers. The assumptions are standard in the literature on fractal processes and are as follows.
Remark 2.2. The technical assumption (A4) is only needed for the asymptotic normality of the estimator of α and not for consistency.
The parameter α ∈ (−1/2, 1/2) is termed the fractal index because it, under mild assumptions, is related to the fractal dimension D = 3 2 − α of the sample paths of the process X (Falconer, 1990;Gneiting et al., 2012). It is also refered to as the roughness index of X, since the value of α reflects itself in the pathwise properties of X, as the following result formalizes.
Proposition 2.1 shows that α controls the degree of (Hölder) continuity of X. In particular, negative values of α corresponds to X having rough paths, while positive values of α corresponds to smooth paths. It is well known that the Brownian motion has α = 0. In Table 1 we give some parametric examples of the kind of processes we have in mind and comment on how they fit into the setup of the present paper; the examples are taken from Table 1 in Gneiting et al. (2012).
To get an intuitive understanding of how the trajectories of the fractal processes look, and in particular how the value of α reflects itself in the roughness of the paths, Figure 1 plots three simulated trajectories of the Matérn process. It is evident how negative values of α correspond to very rough paths, while the paths become smoother as α increases.
The processes in Table 1 are all Gaussian. However, in many applications it is preferable to have rough processes which are both fractal and non-Gaussian (Gneiting et al., 2012, section 3.1.). In the following section we suggest an extension to the above setup that explicitly results in non-Gaussian processes with fractal properties, by considering processes which are volatility modulated.  Table 1, with β = 1, α as indicated above the plots, and n = 500 observations on the unit interval. The same random numbers were used for all three instances.

Extension to stochastic volatility processes
A flexible way to introduce non-Gaussianity of processes for which the theory of the fractal index continues to hold, is through volatility modulation. Following , consider processes of the form where X 0 ∈ R, σ = (σ t ) t≥0 is a stochastic volatility process, and G = (G t ) t≥0 is a zero-mean Gaussian process with stationary increments satisfying (A1)-(A4), e.g. one of the processes from Table  1. The modulation of the increments of G by the stochastic volatility process is a convenient way of introducing non-Gaussianity. To see this, note that the marginal distribution of X t , conditional on the past of the stochastic volatility process and the starting value X 0 , is In other words, the marginal distribution of X t is a normal mean-variance mixture distribution, where the distribution of the stochastic process σ and initial value X 0 determine the mixture. For the integral in (2.2) to be well defined (in a pathwise Riemann-Stieltjes sense), we require that σ has finite q-variation for some q < 1 1/2−α . Intuitively, this means that the "more rough" G is, the "less rough" σ can be. Under these conditions on G and σ, the process X in (2.2) will inherit the fractal properties of the driving process G, as shown in .
For the central limit theorems developed below to hold, we further require another assumption on σ.
(SV) For any q > 0, it holds that for some ξ > 0 and C q > 0.
As pointed out in Bennedsen et al. (2017a), the requirement that σ has finite q-variation for a q < 1 1/2−α can be quite restrictive. For instance, if α < 0 (i.e., G is rough) then σ can not be driven by a standard Brownian motion. A very convenient process, which does not have these restrictions and which is very tractable, is the Brownian semistationary process, which we consider next.

The Brownian semistationary process
Consider X, the (volatility modulated) Brownian semistationary (BSS) process Schmiegel, 2007, 2009), defined as where W is a Brownian motion on R, σ = (σ t ) t∈R a stationary process, and g a Borel measurable function such that t −∞ g(t−s) 2 σ 2 s ds < ∞ a.s. See, e.g., Bennedsen et al. (2017a) for further details of the BSS process. The BSS process is also a normal mean-variance mixture: It is interesting to note that Barndorff-Nielsen et al. (2013) show that for a particular choice of kernel function g and stochastic volatility process σ, X will have a marginal distribution of the ubiquitous Normal Inverse Gaussian type.
We need to impose some technical assumptions on the kernel function g. They are as follows.
(BSS) It holds that (a) g(x) = x α L g (x), where L g is slowly varying at zero.
where L g ′ is slowly varying at zero, and, for any ǫ > 0, we have g ′ ∈ L 2 ((ǫ, ∞)). Also, for some a > 0, |g ′ | is non-increasing on the interval (a, ∞). 2 (c) For any t > 0, The kernel function gives the BSS framework great flexibility. A particularly useful kernel function which has been applied in a number of studies, e.g. Barndorff-Nielsen et al. (2013) and Bennedsen (2017), is the so-called gamma kernel.
Remark 2.3. In Bennedsen et al. (2017b) it was shown that BSS processes satisfying (A1)-(A3), (SV), and (BSS) will have the same fractal and continuity properties as their Gaussian counterparts: for such a BSS process Proposition 2.1 continues to hold. In other words, X will have a modification with Hölder continuous trajectories of order φ for all φ ∈ (0, α + 1/2).

Extension to processes with non-stationary increments
When the increments of X are non-stationary an approach similar to the one in Bennedsen et al. (2017b) can be adopted as follows. Define the time-dependent variogram and, analogously to (2.1), assume that where again C 2,t > 0, α ∈ − 1 2 , 1 2 , and L is a slowly varying function at zero. The methods considered in this paper applies -mutatis mutandis -also to such processes. An example is the truncated Brownian semistationary process.
Example 2.2 (Truncated BSS process, Bennedsen et al. (2017b)). Let where X 0 ∈ R, W is a Brownian motion, and σ a stochastic volatility process. Bennedsen et al. (2017b) call such a process a truncated BSS (T BSS) process. When X satisfies (A1)-(A3) and (BSS), Bennedsen et al. (2017b) show that α is indeed the fractal index of X, in the sense of γ 2 (h, t) satisfying (2.4). We note that processes similar to the T BSS process (with σ t = 1 for all t) have recently been proposed as models of stochastic log-volatility of financial assets, e.g., Gatheral et al. (2018); Bayer et al. (2016).

Summary of assumptions
Above we introduced a number of processes, differing in important ways, most notably through their distributional properties. In spite of these differences, the results presented in this paper will apply equally to all of them. To ease notation, we briefly summarize the assumptions here.
The first set of assumptions is required for consistency of the estimator of the fractal index α.
The second set of assumptions is required for asymptotic normality of the estimator of the fractal index α.
Remark 2.4. As seen from the assumptions, the central limit theorems will not be applicable for α ≥ 1/4. In fact, a central limit theorem do hold in this case, but with a different convergence rate and limiting distribution from what we derive below. When α = 1/4, the convergence rate is n log n and the limiting distribution is zero-mean Gaussian with an asymptotic variance different from when α < 1/4. When α > 1/4 the convergence rate is n 1−2α and the limiting distribution is of the Rosenblatt type, see Taqqu (1979). If one is interested in the range α ∈ [1/4, 1/2) and desire asymptotic normality results similar to what we have below, we recommend using gaps between the observations as in Corcuera et al. (2013) Remark 4.4; the downside of this approach is that one is forced to throw away observations. Given the results presented below, filling in the details of this approach is straight forward, albeit notationally cumbersome. Since the case of very smooth processes, i.e. α ≥ 1/4, seems of limited practical value, we do not pursue this further here.

Semiparametric estimation of, and inference on, the fractal index
Consider n equidistant observations X 1/n , X 2/n , . . . , X 1 of the stochastic process X, observed over a fixed time interval, which we without loss of generality take to be the unit interval, so that the time between observations is 1 n . As n → ∞, this gives rise to the so-called in-fill asymptotics. In what follows, suppose that the process X satisfies the assumptions (A1)-(A3).
When X is Gaussian, it holds, by standard properties on the (absolute) moments of the Gaussian distribution and (2.1), that where p > 0, the function L p (h) := L(h) p/2 is slowly varying at zero, and C p > 0 is a constant. This motivates the regression logγ p (k/n; X) = c p + a log |k/n| + U k,n + ǫ k,n , k = 1, 2, . . . , m, where m ∈ N is a bandwidth parameter, c p = log C p , a = (2α + 1)p 2 , U k,n = log γ p (k/n; X) γ p (k/n; X) , and ǫ k,n = log L p (k/n).
The variogram γ p is estimated straightforwardly aŝ The OLS estimator of the parameter a is naturallŷ with "T" denoting the transpose of a vector and x m being the m × 1 vector while logγ m p := (logγ p (1/n; X), logγ p (2/n; X), . . . , logγ p (m/n; X)) T .
Given an estimateâ OLS of a, our estimate of the fractal index iŝ This estimator is well known and much used in the literature, e.g. Gneiting and Schlather (2004) where " P →" refers to convergence in P-probability.
A number of studies have considered the asymptotic properties of the OLS estimates coming from (3.4), e.g. Constantine and Hall (1994), Davies and Hall (1999), and Coeurjolly (2001Coeurjolly ( , 2008. For a brief summary of this literature, see Gneiting et al. (2012), Section 3.1. The following theorem presents the details in the context of this paper.
Perhaps surprisingly, Theorem 3.1 shows that the asymptotic distribution of the OLS estimator does not depend on the precise structure of the underlying process X, but only on the value of the fractal index α, through the correlation structure of the increments of a fractional Brownian motion (fBm) with Hurst index H = α + 1/2, and possibly the "heteroskedasticity factor" S p . The reason for this is that the small scale behavior of a process X fulfilling assumption (A1), will have the same small scale behavior as increments of the fBm. To see this, write γ 2 ((j + 1)/n; X) − 2γ 2 (j/n; X) + γ 2 ((j − 1)/n; X) 2γ 2 (1/n; X) by assumption (A1) and the properties of slowly varying functions. We recognize (3.7) as the correlation function of the increments of an fBm with Hurst index H = α + 1/2. As shown in the proof of Theorem 3.1, this will imply that the asymptotic variance of the estimator, σ 2 m,p , is the same for all Gaussian processes fulfilling assumptions (A1)-(A4), including the fBm. However, as the theorem also shows, the asymptotic distribution proves to be slightly different when we consider conditionally Gaussian processes. In this case, the stochastic volatility component σ introduces heteroskedasticity, which results in the extra factor S p in the central limit theorem. To make inference feasible in practice, we need to estimate this factor. For this, define whereγ · is given in (3.3). We can prove the following.
where S p is given in (3.6).
Proposition 3.2 shows that S p of (3.8) is a suitable estimator for our purpose: when X is Gaussian, the factor is asymptotically irrelevant, while when X is non-Gaussian (volatility modulated) it provides the correct normalization. This justifies including the factor S p , whether or not one believes the data is Gaussian, at least when any potential non-Gaussianity is volatility induced. In fact, the following corollary is a straightforward consequence of Theorem 3.1, Proposition 3.2, and the properties of stable convergence; the corollary has obvious applications to feasibly conducting inference and making confidence intervals for α.
Corollary 3.1. Suppose the assumptions of Theorem 3.1 hold. Now, where "d" denotes convergence in distribution and σ 2 m,p (α) denotes the asymptotic variance calculated using the estimateα.
Remark 3.2. When using Corollary 3.1 for hypothesis testing, we recommend calculating σ 2 m,p (·) using the value of α under the null, instead ofα.
To apply the above results we need to calculate the factor σ 2 m,p , which boils down to calculating the entries of the matrix Λ p given in equation (3.5). Unfortunately, this is only feasible when p = 2 and becomes increasingly cumbersome as m increases. (The already tedious calculation for p = m = 2 is given in Bennedsen et al., 2016, Appendix B.). For this reason, we recommend Monte Carlo estimation of σ 2 m,p ; in fact, we suggest using the finite sample analogue of this factor. The procedure is detailed in Appendix B; in the next section we present an example of the output, when we study the effect of the choice of bandwidth, m.

Choosing the bandwidth parameter
The choice of bandwidth parameter m is, in general, an open problem. Standard practice in the literature is to set m = 2 (Gneiting et al., 2012, Section 2.3). Indeed, Constantine and Hall (1994) argue that the bias of the estimator increases with m and Davies and Hall (1999) present simulation evidence for the optimal value, in terms of mean squared error, being m = 2. Setting m = 2 amounts to estimating α by drawing a straight line between only the two points closest to the origin, logγ p (1/n; X) and logγ p (2/n; X), when running the OLS regression in (3.2). While tempting from a bias viewpoint, we conjecture that this can result in increased variance of the estimator, by relying on just two points in the regression. In what follows, we examine this in more depth. To be specific, we consider the effect that the bandwidth has on the estimator of the fractal index; first on the theoretical (finite sample) variance ofα, as derived in Theorem 3.1 ( Figure  2), and then on the finite sample bias and mean squared error of the estimator when applied to simulated paths of the various processes of Table 1 (Figure 3). For these investigations, we consider both α = −0.20 (rough case) and α = 0.20 (smooth case). Figure 2 studies the effect that the choice of bandwidth has on the variance of the estimator of α: we plot the approximation of the finite sample variance ofα, σ 2 m,p,n , which is approximately equal to n −1 σ 2 m,p , cf. Theorem 3.1. From the figure, we see that the choice of bandwidth indeed has an effect on the variance of the OLS estimator of α. Interestingly, the effect is very different in the rough case, as compared to the smooth case. In the former, it is evident from the top left plot of Figure 2, that the variance is minimized by an intermediate value of m such as m = 5 or m = 10. To further investigate this, the top right plot shows the ratio between the finite sample variance when m = 2 and when m ∈ {5, 10, 25}. Numbers less than one indicate that the variance of the estimator with m = 2 is greater than the variance of the corresponding estimator with m > 2,  Table  1. We set p = 2 and n = 1 000. The parameter values are given in the text. and vice versa. These ratios seem quite stable as a function of sample size n and it is evident that, from a variance stand point, it is preferable to choose an intermediate m > 2 -indeed, the variance of the estimator is reduced by approximately 40% when going from m = 2 to m = 5. These conclusions get turned on their heads when we consider the smooth case, α = 0.20, in the bottom row: here it seems that m = 2 is optimal. We further investigate this through simulations as in Davies and Hall (1999): Figure 3 plots the bias (left) and mean squared error (right) of the estimator (3.4), as a function of bandwidth m, for the five parametric processes of Table 1. To calculate the finite sample bias and mean squared error of the estimator, we simulate B = 10 000 instances of each process, each with n = 1 000 observations; the true value of the fractal index in this exercise is α = −0.2 (top row) and α = 0.2 (bottom row). The scale parameter is set to β = 1. For the Cauchy and Dagum processes we additionally set τ = 1 and τ = 0, respectively. When looking at the rough case, α = −0.2, the same conclusion as above emerges: even though the bias do increase, as expected, with increasing m, it is clear that the mean squared error is minimized for an m > 2. In this case, i.e. for these parameter values and this sample size, the minimum is attained between m = 5 and m = 10 for all five processes. We again conclude that an intermediate value for the bandwidth is preferable in finite samples when α < 0. The smooth case, α = 0.2, also matches what we found above: indeed, we find that both bias and mean squared error increase with increasing m, so here m = 2 seems optimal.
In conclusion, the evidence of this section suggests that when the underlying process is rough, the optimal choice of bandwidth is some m > 2 and we recommend an intermediate value such as m = 5. In contrast, when the process is smooth, m = 2 is preferable. Although setting m = 2 seems to be accepted practice in the literature, we believe that the rough case of α < 0 is arguably more relevant in empirical applications. For this reason we suggest using an intermediate value for the bandwidth parameter, unless one has reason to believe the underlying data to be smooth. 3

Asymptotic theory in the presence of additive noise
Consider now the situation, where the observations of X, satisfying (A1)-(A3), are contaminated by additive noise; that is, instead of observing X, we observe the process Z, given by Z j/n := µ + X j/n + u j , j = 1, 2, . . . , n, ( 3.9) where µ ∈ R is a constant and u = {u j } n j=1 is a Gaussian iid noise sequence with mean zero and variance σ 2 u := V ar(u 1 ) ≥ 0. (When σ 2 u = 0 we mean that the noise is absent from the observations) Since we observe Z, and not X, what is relevant for us is the "contaminated", or "noisy", variogram, i.e. the variogram of the observation process Z: where the last equality follows from Assumption (A1). From this we see that when σ 2 u > 0, log γ 2 (h; Z) will not be linear in log h, hence the estimator (3.4) of α will not be applicable; in fact, it is not hard to show that this estimator will be downwards biased in the presence of noise, i.e. when applied to γ p (·; Z). In fact, the following is true.
Proposition 3.3. Suppose that the observations of a process Z are given by (3.9) with σ 2 u > 0, where X satisfies assumption (LLN). Fix p > 0, m ∈ N, and letα =α p,m be the OLS estimator of α from (3.4) using the contaminated version of the empirical variogramγ p (·; Z) in place ofγ p (·; X) in the regression (3.2). Now, Proposition 3.3 shows that if the data are contaminated by noise, then estimates of the parameter α will be biased downwards towards −1/2, i.e. the lowest permissible value for α. In other words, if the data are contaminated by noise, then the estimator of α considered above, will lead one to conclude that the data are more rough than what is actually the case for the underlying process X. This is an important point to note for the practitioner: when finding evidence of roughness (i.e. α < 0) in data, it is crucial to consider whether this is due to an intrinsic property of the underlying data generating mechanism or whether it could simply be the product of noise, e.g. measurement noise.
Fortunately, it is possible to account for the noise when estimating α to arrive at a consistent estimator. For instance, Bennedsen et al. (2017a) suggest a noise-robust estimator based on a non-linear least squares regression -however, this estimator does not allow for the slowly varying function L and requires the interval over which the process is observed to grow. Presently, therefore, we propose an alternative noise-robust estimator which is valid in our in-fill asymptotics setup and again relies on a simple OLS regression.
First, for an integer κ ≥ 2, define the function From (3.10) and assumption (A1), we have where it is easy to show that the function is slowly varying at zero. From this, it is clear that the logarithm of f p (h; Z, κ) is -up to the slowly varying function L * p -linear in log h. This motivates a linear regression as the one in (3.2) with logf p in place of logγ p : logf p (k/n; Z, κ) = b * + a * log |k/n| + U * k,n + ǫ * k,n , k = 1, 2, . . . , m, (3.12) wheref p (k/n; Z, κ) :=γ p (κk/n; Z) 2/p −γ p (k/n; Z) 2/p is the empirical estimate of the function f , which is feasible to calculate from the observations Z j/n . Define the noise robust estimate of α aŝ whereâ * OLS is the OLS estimate of a * from the linear regression (3.12), analogous to (3.4) withf p in place ofγ p . We can prove the following.
Remark 3.3. Proposition 3.4 allows for σ 2 u = 0 i.e. for there to be no noise in the observations. In other words, the robust estimator is a consistent estimator of α, also in the absence of noise. In Figure 4 we illustrate the use of Proposition 3.4 by calculating the bias and root mean squared error (RMSE) of the two OLS estimators given in (3.4) and (3.13), when applied to a process Z with σ 2 u > 0. The details are provided in the caption of the figure. The former estimator is not robust to the noise in Z, while the latter estimator is per Proposition 3.4. It is clear how this manifests itself in a large bias in the OLS estimator (3.4). In fact, although the true value of the fractal index of the underlying process is α = −0.20, the mean OLS estimates coming from the non-robust estimator is −0.4608, i.e. almost at the lowest permissible value of −1/2. This is of course a consequence of Proposition 3.3. In contrast, the robust estimator (3.13) proposed in this section is practically unbiased for most values of the parameter κ, at least when κ ≥ 4.
Although the results of this section hold for all integer κ ≥ 2, the actual finite sample performance of the results can be quite sensitive to this tuning parameter, as also witnessed in Figure  4. The optimal choice of κ seems to depend on the number of observations n and the variance of the noise σ 2 u ; an investigation into the exact way this is the case is beyond the scope of the present paper. In practice, we recommend that the researcher run some numerical experiments on simulated data under conditions similar to those of the practical experiment; simulation experiments such as the one in Figure 4 for example. We provide an example of how one can construct such a simulation experiment to arrive at a reasonable value for κ in Section 5.2, where we apply the robust estimator to a time series of financial prices.
The next result provides the central limit theorem, as it relates to the robust estimator.
where S p is as in Theorem 3.1 and σ 2, * m,p is given in Appendix C.
Remark 3.4. A feasible central limit theorem in the case of (i) is straightforwardly constructed in the same way as in Corollary 3.1, including Monte Carlo estimation of σ 2, * m,p . As shown in (ii) of Theorem 3.2, the presence of the noise will unfortunately result in a variance ofα * , which decays slower than √ n; indeed, the exact distribution ofα * is difficult to derive and even harder to feasibly estimate.

A test for the presence of noise
Using the above, we can now construct a test for whether the observed time series Z contains noise or not. To be specific, we are interested in testing the null hypothesis H 0 : σ u = 0 against the alternative H 1 : σ u > 0. (3.14) Tests of this kind, in the context of time series of asset prices, were considered in Aït-Sahalia and Xiu (2018), where the authors develop a test for the presence of market microstructure noise in high frequency data. The test proposed here is similar in spirit to the test of Aït-Sahalia and Xiu (2018) and in Section 5.2 we briefly consider testing for the presence of market microstructure in high frequency asset prices as well.
To device the test, we consider the difference between the robust estimatorα * from (3.13) and the usual (non-robust) estimatorα from (3.4). From Propositions 3.1 and 3.4 is is immediately clear that under H 0 while under H 1 , Proposition 3.3 additionally implies that Analogously to Theorem 3.2, we can also prove the following.
where S p is as in Theorem 3.1 and σ 2, * * m,p is given in Appendix C. . (3.15) The following corollary is a straightforward application of Theorem 3.3.
Corollary 3.2. LetÂ n be as in (3.15). Now, The applicability of Corollary 3.2 for testing whether a fractal process is contaminated by noise is obvious.
Remark 3.5. Above we have assumed that the noise sequence u is Gaussian. However, one can show that all the results of sections 3.2 and 3.2.1 apply for general iid noise sequences with finite variance when p = 2. In other words, if the Gaussian assumption on the noise sequence u is not fulfilled -or seems too restrictive -then one should choose p = 2 and go ahead and apply the results of these sections.

Simulation studies
To examine the finite sample properties of the central limit results presented above, we here conduct three small simulation studies and collect the results in Tables 2-4. In each study we will let X be an fBm with Hurst index H = α + 1/2, for various values of α, and simulate n observations on the interval [0, 1]. Additional information on the exacts simulation setups are given in the captions of the tables.
For a value α 0 ∈ (−1/2, 1/4), Corollary 3.1 allows us to test the null hypothesis Panel A of Table 2 shows the empirical size -i.e., the rejection rates of H 0 when H 0 is true -of this test for various values of n and α, at a nominal 5% level. In contrast, Panel B shows the empirical local power of the test, that is, the rejection rate of the null hypothesis H 0,n : α = α 0 + n −1/2 against the alternative H 1,n : α = α 0 + n −1/2 , (4.2) when α 0 is the true value of α used in the simulations of the fBm and n is the number of observations. In other words, the value of α that we test for asymptotes towards the true value of α at the rate n 1/2 . We see that both the empirical size and power properties of the test are quite good. In particular, the test is approximately correctly sized for n ≥ 20. (We conjecture that some of the already small deviation from the nominal 5% level is likely due to Monte Carlo simulation error.) Similarly, the empirical power is also adequate. We conclude that the test is very precise and can comfortably be used in practice, even for small sample sizes. Table 3 contains the analogous results in the setup of Theorem 3.2 (i). Note that we here set σ 2 u = 0 so that we do not have any noise in the observations. We set κ = 10 for all sample sizes but this could likely be optimized by considering different values of κ for each different value of n. In this way, the empirical properties of the test presented in Table 3 are conservative, in the sense that one could likely obtain even better properties by choosing κ in a data-driven way, cf., e.g., the approach taken in the setup of Figure 4. When considering the local power in this case, we replace the rate n 1/2 by n 1/4 , that is, in Panel B of Table 3 we are testing the null hypothesis H 0,n : α = α 0 + n −1/4 against the alternative H 1,n : α = α 0 + n −1/4 , (4.3) where α 0 is the true value of α used in the simulations of the fBm. Again we conclude that the empirical properties of this test are good. We see that we need more observations here to get a properly sized test, than what was the case above, the reason being that the estimator now utilizes gaps of size κ = 10 between observations, when calculating the statisticf . The empirical power properties also seem adequate, although the rate is now somewhat slower than what we saw above. All in all, though, the practical relevance of this hypothesis test is not as great as the one considered above: we do not posses the asymptotic theory in the case noisy case (σ 2 u > 0) and in the non-nosiy case (σ 2 u = 0) the setup of the regular OLS estimator is preferable for inference. The robust estimator is consistent, though, so for pointwise estimation, the estimator should still have great practical relevance.
Finally, Table 4 contains the results from the test (3.14); that is, we test for the presence of noise in the observation using Corollary 3.2. Now the true observations are Z i = µ + X i + u i , where the variance of u i is σ 2 u ≥ 0. In Panel A we again have the size of the test, which is where we set σ 2 u = 0, i.e. where H 0 is true. In Panel B we have the power properties, which is the rejection rates of H 0 when σ 2 u = 0.05, i.e. when H 1 is true. 4 We see that the test is slightly under-sized in most cases, although reasonably close to the nominal level. The power properties are reasonably good for n ≥ 800 at least when α is not too negative; when α ≈ −0.50, the test apparently needs many observations to be able to reject the null when it is false. This is to be expected, however, considering Proposition 3.3. Overall, the empirical properties of the test are quite good and we remind the reader that they can be further improved by choosing κ in a data-driven way, see 4 As can be seen in Table 4, in this simulation study we consider fewer Monte Carlo replications than what we did in the two studies above; this is because we here, contrary to above, do not have a null value for α -hence we need to use the estimated value,α * , when calculating the asymptotic variance by Monte Carlo simulation (cf., Appendix B), whereas we in the other two studies could use the null value when calculation the asymptotic variance of the estimator. The upshot is, that in the present study, we have a Monte Carlo simulation (of the asymptotic variance given the estimateα * ) inside a Monte Carlo simulation (of the empirical rejection rates of the hyptohesis test), making the computational burden rapidly increasing in the number of Monte Carlo replications. We therefore choose a lower value for this number in the third study. In practice, for a given time series, one only needs to run one Monte Carlo study (approximating the asymptotic variance given the estimateα * ) which is quite fast -hence in practice, the test is still feasible to do at high accuracy, i.e. with many Monte Carlo replications in the estimation of the asymptotic variance.

Application to turbulent velocity data
We first illustrate the use of the above methods in the context of the study of turbulent velocity flows. More precisely, we have at our disposal a time series of one-dimensional measurements of the longitudinal component of a turbulent velocity field in the atmospheric boundary layer, measured 35 meters above ground level. The series consists of 20 million equidistant observations over a time period of 4 000 seconds, i.e. with 5 000 measurements per second (in other words, a sampling frequency of 5 kHz). The same time series has been studied in e.g. Corcuera et al. (2013), Barndorff-Nielsen et al. (2014), Bennedsen et al. (2016), and we refer to Dhruva (2000) for more information on the data set itself. We first de-mean the long time series and standardize it to have unit variance. We then sample the data at 5 Hz, i.e. 5 observations per second; this frequency is squarely in the so-called inertial range, where the celebrated 5/3-law of Kolmogorov (1941a,b) states that -in the context of this paper -the process underlying the time series should have fractal index α = −1/6 (Corcuera et al.,

2013, Section 5).
The left plot of Figure 5 contains the standardized data series sampled at 5 Hz. The right plot gives examples of the regression (3.2) for p ∈ {1, 2, 3}. In this plot, the crosses are logγ p (k/n; X) as a function of log(k/n), while the lines are the associated OLS regression fits using a bandwidth of m = 5. The estimated slopes of the three lines are respectivelyâ OLS = 0.3361, 0.6713, and 1.0044 which, through (3.4), yield the following estimates of the fractal index,α = −0.1643 (0.0050), −0.1638 (0.0047), and −0.1630 (0.0049), where the numbers in parentheses denote the standard deviation of the estimates, cf. Theorem 3.1. Likewise, using Corollary 3.1, we test the the prediction of Kolmogorov via the null hypothesis H 0 : α = −1/6 against the alternative H 1 : α = −1/6, at a nominal 5% level. The corollary in this case yields P-values of 63.38%, 54.52%, and 45.54%, respectively. In other words, we can not reject the prediction that α = −1/6 in this data set of turbulent velocity flows at this sampling frequency.

Application to financial price data
An oft-used model of high frequency financial logarithmic prices is Z j/n = X j/n + u j , j = 1, . . . , n, where u = {u j } n j=1 is a market microstructure noise process (O'Hara, 1995) and X is a stochastic volatility process, e.g., It is well known that for the absence of arbitrage in the market, it is a requirement that X is a semimartingale (Delbaen and Schachermayer, 1994), which under mild assumptions, is equivalent to G, and therefore X, having α = 0. In other words, since we expect the market to be free of arbitrage, we would expect that X has fractal index α = 0.
To test this, we study a long time series of the logarithm of financial futures prices, supposing they come from the model (5.1). To be specific, we have at our disposal data recorded every second on the front month E-mini S&P 500 futures contract, traded on the CME Globex electronic trading platform, from January 3, 2005 until December 31, 2014. We exclude weekends and holidays and keep only full trading days, which results in 2 495 days. For these days, we further restrict attention to the most active period of the day, which is when the New York Stock Exchange (NYSE) is open, from 9.30 a.m. until 4 p.m. Eastern Standard Time (EST). This results in 23 400 seconds (6.5 hours) for each day. We estimate α each day using both estimators (3.4) and (3.13); this results in N = 2 495 estimates of α each calculated from n = 23 400 observations.
As discussed, since we believe the market to be free of arbitrage opportunities, we would expect to find α ≈ 0 for financial price series; however, when we estimate α from the data using the standard OLS estimator (3.4) we most often find very negative values; in fact the mean estimate of α over the N = 2 495 days is −0.20. Similarly, when we use Theorem 3.1 to test the null hypothesis H 0 : α = 0 against the alternative H 1 : α = 0, (5.2) at a nominal 5% level, we reject H 0 on 98.92% of the days. In other words, at first glance, it seems that high frequency log prices are very rough. Of course, this finding might simply an artifact of the noise sequence u which we expect to be present in high frequency stock prices; we therefore also apply the robust estimator (3.13). The robust estimator requires a choice of tuning parameter κ and, as suggested in Section 3.2, we run a simulation experiment to gauge a reasonable value for this parameter. We set up our experiment to be realistic so that we can expect the optimal value of κ we find to also be a good value to use on the real data. We therefore simulate B = 5 000 instances ("days") of an fBm with n = 23 400 observations ("seconds"); for these fBms we set α = 0, as we expect this to be the true value from the underlying process. (when α = 0 the fBm is a standard Brownian motion.) For each simulated fBm, we add an iid sequence of Gaussian noise with σ 2 u = 0.01; this value of σ 2 u was chosen since it is a reasonable value for the variance of microstructure noise (Hansen and Lunde, 2006). 5 We then apply (3.13) and calculate the RMSE for various values of κ ≥ 2. The results are shown in the left plot of Figure 6 where we see that the RMSE is minimized for κ ≈ 60.
We then go ahead and apply both estimators to the empirical data described above; the results are seen in the right plot of Figure 6. As just described, the non-robust OLS estimator (blue dots) are often very negative (the average across all days is −0.20). Conversely, the robust estimator with κ = 60 (red crosses) seem much closer to the expected value of zero (the average across all days is −0.0448). We also apply Corollary 3.2 to the data, to test the null hypothesis, (3.14), of no noise in the observations. When performing this test each day at a 5% level, we reject the null on 99.72% of the days. In other words, the formal test provides very strong evidence for the presence of market microstructure noise in the observations.

Conclusion and open problems
In this paper, we have laid out a large and coherent framework for analyzing data with fractal characteristics. We focused on a particular estimator of the fractal index but the methods can be applied more generally and in particular to other related estimators as well. For instance, the extension to using higher order differences -instead of first order differences, as we used herewhen calculating the variogram, is straight forward. We showed how to extend the theory to a large class of non-Gaussian processes. The consistency result of the estimator turns out to be the same for Gaussian and non-Gaussian processes, while the central limit theorem requires an extra factor in the non-Gaussian case. The estimator of this correcting factor -S p in Corollary 3.1 -is asymptotically irrelevant when no stochastic volatility is present. For this reason we recommend always including this factor when applying the central limit theorem of the estimator of the fractal index, whether or not one believes the underlying data to be Gaussian. Lastly, we saw that noise in the observations, e.g. measurement noise, will bias the estimates of the fractal index downwards. It is important that the practioner is cognizant of this possible bias when studying data which are potentially rough; we suggested an estimator which is robust to such noise and can be applied if one expects the data to be so contaminated.
Let us briefly comment on a few possible directions for further study. The OLS estimator relies on a bandwidth parameter m; in Davies and Hall (1999) the authors find, through simulations, that the optimal value for this is m = 2. However, as we saw above, this was not the case in general when we studied the finite sample properties of the estimator. Indeed, for rough processes, an intermediate value for the bandwidth, such as m = 5, is preferable. How to choose m optimally in finite samples is an open and interesting problem. Likewise, what value of the power parameter p > 0 to use is an open problem; Gneiting et al. (2012) recommend p = 1 as they find that this makes the estimation more robust to, e.g., outliers or non-Gaussianity. However, above we proposed an estimator which is robust to such anomalies by construction, in particular when p = 2. An in-depth investigation of the properties of the robust estimator as compared to the standard estimator with p = 1 would valuable. Lastly, we also saw that the robust estimator relies on the tuning parameter κ; finding a data-driven way to choose this parameter would be useful. In general, the problem of "robust estimation" is interesting and more theoretical as well as applied research into this field would be valuable.
To see that the term x T m ǫ m vanishes as n → ∞, note that and therefore x T m ǫ m = m k=1 x m,k log L p (k/n) = m k=1 x m,k log L p (k/n) L p (1/n) → 0, n → 0, since lim n→∞ Lp(k/n) Lp(1/n) = 1 by the property of slowly varying functions. The required result now follows by noting that γ p (k/n; X) γ p (k/n; X) =γ p (k/n; X) m p γ 2 (k/n; X) p/2 P → 1, n → ∞, k ≥ 1, which holds by Propositon 1 in   where Λ p = {λ k,v p } m k,v=1 is a m × m matrix with entries λ k,v p = lim n→∞ n · Cov γ p (k/n; B H ) γ p (k/n; B H ) ,γ p (v/n; B H ) γ p (v/n; B H ) , k, v = 1, 2, . . . , m, (A.4) with γ p (·; B H ) denoting the p'th order variogram for a fractional Brownian motion with Hurst index H = α + 1/2, and similarly forγ p . Note that the limit in (A.4) exists for k, v = 1, 2, . . . m, by Breuer and Major (1983), Theorem 1, see also Corcuera et al. (2013), Remark 3.3. As we will show below, assumption (A4) implies that √ nx T m ǫ m → 0, n → ∞, (A.5) which means that from (A.1) we get, using (A.3) and the delta method, which is what we wanted to show. To see that (A.5) holds, use the rule of l'Hôpital, and the properties of slowly varying functions, to conclude lim n→∞ √ n log L p (k/n) L p (1/n) = lim n→∞ 2n −1/2 L ′ p (k/n) L p (k/n) k − L ′ p (1/n) L p (1/n) = 0, by assumption (A4). This concludes the proof when X is Gaussian. Suppose instead that assumption (CLTc) holds (the case (CLTb) is similar). We proceed analogously to above: by Theorem 4 of Barndorff-Nielsen et al. (2011), see also Theorem 3.2. and Remark 3.4. of Corcuera et al. (2013), we get where B is an m-dimensional Brownian motion, defined on an extension of the original probability space, (Ω, F, P), independent of F. The matrix Λ p is identical to the one above. We proceed as before. In particular, invoking the delta method we get and where Z p is as in Theorem 3.1. This concludes the proof.
The result now follows from Proposition 1 in .
Proof of Proposition 3.3. This follows easily from the fact that when σ 2 u > 0 it holds for any k ≥ 1 γ p (k/n; Z) P → γ p (0; Z) = C p σ p u > 0, n → ∞, so that

A.1 Proofs related to the noise robust estimator
Most of the proofs related to the noise robust estimatorα * proceeds analogously to the proofs related to the standard estimatorα given above. Indeed, the function f p (h; Z, κ) has the same behavior as γ 2 (h; X) for h ≈ 0 -cf. equations (3.1) and (3.11). The upshot is that most of the theory will go through as in the proofs above. The following Lemma formally provides the details.
With this function L * (2) , it is also the case that (A3) holds. To see this, write L * L(x) .
Using that both L and L 2 are slowly varying at zero and obey (A3) by assumption, the result follows.
B Monte Carlo approximation of the variance of the OLS estimator