A Generalized Correlated Random Walk Converging to Fractional Brownian Motion

TL;DR

This paper introduces a new algorithm for generating fractional Brownian motion with a specified Hurst parameter using correlated Bernoulli variables, establishing its convergence to the process through theoretical proof.

Contribution

The authors develop a novel correlated Bernoulli-based algorithm for fractional Brownian motion generation and prove its convergence to the target Gaussian process.

Findings

The proposed method accurately generates fractional Brownian motion with specified Hurst parameters.

Theoretical proof confirms the convergence of the random walk to fractional Brownian motion.

The algorithm leverages the correlation structure between Gaussian and Bernoulli variables.

Abstract

We propose a new algorithm to generate a fractional Brownian motion, with a given Hurst parameter, 1/2<H<1 using the correlated Bernoulli random variables with parameter p; having a certain density. This density is constructed using the link between the correlation of multivariate Gaussian random variables and the correlation of their dichotomized binary variables and the relation between the correlation coefficient and the persistence parameter. We prove that the normalized sum of trajectories of this proposed random walk yields a Gaussian process whose scaling limit is the desired fractional Brownian motion with the given Hurst parameter, 1/2<H<1

Equations

$$E[B^H(t)B^H(s)]=\frac{1}{2}(t^{2H}+s^{2H}-|t-s{|}^{2H}).$$

$$h^H(m)=\frac{1}{2}[(m+1{)}^{2H}+(m-1{)}^{2H}-2m^{2H}],$$

$${\sigma }_{jk}=\frac{\Phi [z(p_j),z(p_k),{\delta }_{jk}]-p_j{p}_k}{\sqrt{p_j(1-p_j)p_k(1-p_k)}},$$

$$\Phi [z(p_j),z(p_k),{\delta }_{jk}]={\int }_{-\infty }^{z(p_j)}{\int }_{-\infty }^{z(p_k)}f(z(p_j),z(p_k),{\delta }_{jk})dz(p_j)dz(p_k)$$

$${\sigma }_{jk}=\frac{\Phi [z(p),z(p),{\delta }_{jk}]-p^2}{p(1-p)},$$

$$\Phi [z(p),z(p),{\delta }_{jk}]={\int }_{-\infty }^{z(p)}{\int }_{-\infty }^{z(p)}f(z(p),z(p),{\delta }_{jk})dz(p)dz(p)$$

$${\delta }_m^H=\frac{1}{2}[(m+1{)}^{2H}+(m-1{)}^{2H}-2m^{2H}],$$

$${\sigma }_{n,n+1}=\frac{\Phi [z(p),z(p),{\delta }_1^H]-p^2}{p(1-p)},$$

$${\sigma }_{n,n+m}=\frac{(4\Phi [z(p),z(p),{\delta }_1^H]-4p+1{)}^m-(2p-1{)}^2}{4p(1-p)}$$

$$E[{ϵ}_m^{\rho }{ϵ}_{m+n}^{\rho }]=(2\rho -1{)}^n$$

$$E[{ϵ}_m^{\mu }{ϵ}_{m+n}^{\mu }]={\int }_0^1(2\rho -1{)}^{n}d\mu (\rho )$$

$$\rho =2\Phi [z(p),z(p),{\delta }_1^H]-2p+1,$$

$$P(Y_{n+1}^{\rho }=Y_n^{\rho }|\sigma (Y_k^{\rho },0\le k\le n))$$

$$P[Z_n\le z(p),Z_{n+1}\le z(p)]=\Phi [z(p),z(p),{\delta }_{n,n+1}].$$

$$P[Z_n>z(p),Z_{n+1}\le z(p)]$$

$$P[Z_n>z(p),Z_{n+1}>z(p)]$$

$$\rho =2\Phi [z(p),z(p),{\delta }_1^H]-2p+1$$

$$r(Y_m,Y_{m+n})=\frac{(4\Phi [z(p),z(p),{\delta }_1^H]-4p+1{)}^n-(2p-1{)}^2}{4p(1-p)}$$

$$r(Y_n,Y_{n+1})=\frac{\Phi [z(p),z(p),{\delta }_1^H]-p^2}{p(1-p)}=\frac{E[Y_n{Y}_{n+1}]-(2p-1{)}^2}{4p(1-p)}$$

$$E[Y_{n+1}|\sigma (Y_k,0\le k\le n)]$$

$$r(Y_m,Y_{m+n})$$

$$g(p)$$

$$\underset{N\to \infty }{\lim }\underset{M\to \infty }{\lim }{a}_H\frac{X_{[Nt]}^{H,1}+X_{[Nt]}^{H,2}...+X_{[Nt]}^{H,M}-M(Nt)(2p-1)}{N^H\sqrt{M4p(1-p)}}=B^H(t)$$

$$f(\rho )=(1-H)2^{3-2H}(1-\rho {)}^{1-2H}$$

$$E((X_k^{H,i}{)}^2)$$

$$\underset{M\to \infty }{\lim }\frac{X_k^{H,1}+X_k^{H,2}...+X_k^{H,M}-Mk(2p-1)}{\sqrt{M}\sqrt{4p(1-p)}}=Z_k^H$$

$$G_k^H:$$

$$r(n)$$

$${\sigma }_{k,k+n}$$

$$g(p)$$

Full text

Keywords Correlated Random Walk, Dichotomized Binary Variables, Fractional Brownian motion, Gaussian Process

1 Introduction

Most of the real data displaying long-range dependence can be modeled using self-similar processes. Fractional Brownian motion (fBm) is one of the simplest models demonstrating long-range dependence. That is why this phenomenon has been widely used in many research areas and disciplines in recent years. For instance, in communication systems, Leland et al., (1994) used the increments of fBm for modeling Ethernet local area network (LAN) traffic. In the field of mathematical finance, Rogers, (1997) proposed an fBm model to explain the movement of share prices. In biology, Lim and Muniandy, (2001) utilized the discrete-time version of fBm to model the non-coding sequence of human DNA by recognizing a DNA sequence as a fractal random walk. Additionally, there is still an increasing interest in the use of fBm in many fields since it well captures the dependence behavior of the data. For this reason, its theoretical properties and path behaviors have been well studied and the interest in these studies still continues.

Definition 1

Let $H$ be a constant belonging to $(0,1)$. A fractional Brownian motion (fBm), $B^{H}(t)_{t\geq 0},$ is continuous, centered Gaussian process with Hurst parameter $H,$ and with covariance function;

[TABLE]

The fractional Brownian motion satisfies the following properties;

1. i.

$B^{H}(0)=0$, and $E(B^{H}(t))=0$ for all $t\geq 0.$ 2. ii.

$B^{H}(t)$ has stationary increments that is $B^{H}(t+s)-B^{H}(s)$ has the same distribution with $B^{H}(t)$ for $s,t\geq 0$. 3. iii.

$B^{H}(t)$ is a Gaussian process. 4. iv.

The variance of $B^{H}(t)$ equals $t^{2H}$ for all $t\geq 0$ and $H\in(0,1)$. 5. v.

$B^{H}(t)$ has continuous trajectories. 6. vi.

The increment process of fBm, $\{B^{H}(n+1)-B^{H}(n):n=0,1,2,...\},$ also called the fractional Gaussian noise(fGn) which are jointly Gaussian variates with zero mean and autocovariance function

[TABLE]

where $H$ is the Hurst index describing the dependency among the increments and integer valued $m,$ Biagini et al., (2008).

It is already known that these increments can be either positively or negatively correlated depending on the value of the Hurst parameter, $0<H<1$. In particular, an fBm with parameter $H=1/2$ corresponds to a standard Brownian motion which has independent increments. For $H<1/2$, its increments are negatively correlated and they display short-range dependence. In contrast, for $H>1/2$, the auto-covariance of the fBm increments is positive. Thus, the two consecutive increments tend to have the same directions, Biagini et al., (2008).

The simulation of fBm is of crucial importance for its application in fields including economics, finance, engineering and hydrology. Therefore, the development of an algorithm to simulate an fBm is both theoretically and practically required. In the literature, there is a large number of simulation methods. For instance, the integral representation introduced by Mandelbrot and Van Ness, (1968) is used for a direct approximation of fBm. In another study, Hosking, (1984) generates an fBm by developing the method that implicitly computes the fGn covariance matrix. Another approach is the fast Fourier transform (FFT) method developed by Davies and Harte, (1987). In this method, FFT algorithm is used in order to generate an fGn sample. Then, the covariance matrix of fGn is embedded in a circulant covariance matrix, and this circulant matrix is diagonalized with an FFT algorithm. The Cholesky method proposed by Asmussen, (1998) relies on the Cholesky decomposition. The other approximate and fast technique is the random midpoint displacement method proposed by Lau et al., (1995). Similar to the Hosking, (1984) method, it utilizes an approach based on computing the conditional distribution of fGn. The only difference is that the generation is based on the conditional distribution given certain points instead of all points. Another way of generating an fBm process is the wavelet based synthesis method, see Wornell, (1990), Abry and Sellan, (1996). The principle idea of this method is to write an fBm as a weighted sum of orthonormal wavelets, where the weights are samples of independent centered Gaussian processes. The fast and efficient implementation of this method is given in Abry and Sellan, (1996). The method firstly generates the wavelet coefficients corresponding to an orthogonal basis. Afterwards, fBm is obtained via an inverse wavelet transformation. Caglar, (2000) summarizes and discusses these techniques and introduces a new algorithm using micropulses approximation to synthesize a fractional Brownian motion.

In this study, different from those reviewed above, we propose a new fBm simulation method through a discrete process, namely correlated random walk. To the best of our knowledge, the generation methods of fBm using the discrete processes go back to Dasgupta, (1997), which uses the independent binary random variables and the stochastic integral representation of fBm to approximate an fBm. On the other hand, Sottinen, (2001) defines a random walk, for long range dependence case which converges weakly to fBm using a kernel function that converts the standard Brownian motion to fBm. Szabados, (2001) benefits from the moving average of an appropriately nested sequence of random walks uniformly converging to fBm when $H\in(\frac{1}{4},1)$. This approximation uses the discrete form of moving average representation. As another way of fBm generation by discrete processes, Enriquez, (2004) proves that the normalized persistent random walk converges weakly to fBm. The construction relies on the correlated jumps in such a way that the probability of having the same jump as the previous one, which is the persistence parameter, defines the correlation of random walk. Konstantopoulos and Sakhanenko, (2004) have introduced that scaled random walk using the weighted sum of independent and identically distributed random variables converge to fBm under the sufficient condition for the weak convergence of normalized sums to fBm with $H>\frac{1}{2}$. Similar to Konstantopoulos and Sakhanenko, (2004), Lindstrøm, (2007) uses the same approximation for the case $H<\frac{1}{2}$.

In the Donsker, (1951) Theorem it is proven that standard Brownian motion can be constructed by independent increment random walk using the Central Limit Theorem. As an analogue of this idea, fBm can also be constructed by dependent increment random walk. Except for Enriquez, (2004) the above mentioned studies use independent increment random walk. Enriquez, (2004) construction depends on the persistent random walk with a persistence parameter corresponding to the probability of producing the same jump with the previous one.

In this study, we propose a new algorithm, where we generate correlated binary random variates with proportions $p.$ An explicit density is assigned to the values of $p$ using the link between the correlation of multivariate Gaussian and the correlation of their dichotomization along with the relation between the correlation and the persistence parameter given in Enriquez, (2004). The correlations of increments of this random walk are obtained from the one correlation of increments of fBm for a given Hurst parameter and the marginal proportions, $p.$ We prove that the normalized sum of independent trajectories of this newly proposed correlated random walk yields a discrete Gaussian process by the Central Limit Theorem. Furthermore, its scaling limit is the desired fractional Brownian motion with $H\in[1/2,1],$ since it owns the correlation structure which satisfies the conditions specified in Taqqu, (1975). Our newly proposed algorithm generalizes the construction given in Enriquez, (2004), Thm 1, in a way that it uses the same correlation of increments, but with the Bernoulli $p$ random variables as increments where $p$ has a certain density function, rather than using only Bernoulli $p=1/2$ random variables. This newly proposed algorithm is simple to follow and implement by a broad range of researchers working in various fields of applied sciences and it also has theoretical interest as it is the exact analogous of independent random walk converging to Brownian motion for dependent random walk converging to fractional Brownian motion.

2 The correlated random walk construction

Different from binary variables obtained by natural means such as male or female, yes or no, success or failure, some dichotomous variables can be generated via discretization of continuous ones. The dichotomous variables, such as having a high or low income or being short or tall, can be produced by assigning a threshold value to continuous variables. Despite the loss of information, these variables are crucial for many scientific fields and they are easy to implement and use in many areas such as psychology, criminology, biology and sociology. The correlation between the two continuous variables is usually calculated by Pearson correlation, but when the two are dicotomized by a threshold term, the correlation name changes. The tetrachoric correlation coefficient is assigned for the correlation between two dichotomous variable before discretization and after discretization. The correlation between the dichotomized variables is called phi correlation coefficient. Demirtas and Vardar-Acar, (2017) emphasize that when both variables are discretized, the magnitude of these correlations can easily be transformed into the binary case under the normality assumption and the connection between the tetrachoric correlation and the phi coefficient is known.

Assume that $Z_{j}$’s denote the Normally distributed variables. These are dichotomized to produce $Y_{j}$’s which represent binary variables. Suppose $Y_{1},Y_{2},...,Y_{j}$ are $J$ binary random variables having the mean $E[Y_{j}]=p_{j}$ for $j=1,2,....,J$ and the correlation $Corr[Y_{j},Y_{k}]=\sigma_{jk}$ for $j=1,2,...,J-1;k=2,3...,J$. Let $Z=(Z_{1},Z_{2},...,Z_{J})^{T}$ represent the $J$-dimensional multivariate normal random variables with zero mean and $Corr[Z_{j},Z_{k}]=\delta_{jk}$ for $j=1,2,...,J-1;k=2,3...,J$. Then, the link between the tetrachoric correlation ($\delta_{jk}$) and the phi coefficient ($\sigma_{jk}$) is as follows:

[TABLE]

where $z(p)$ represents $p$ th quantile of standard normal distribution and

$\Phi[z(p_{j}),z(p_{k}),\delta_{jk}]$ be a cumulative distribution function of standard bivariate normal with correlation coefficient $\delta_{jk}$ for $j=1,2,...,J-1;k=2,3...,J.$ Explicitly,

[TABLE]

where $f(z(p_{j}),z(p_{k}),\delta_{jk})=[2\pi^{-1}(1-\delta_{jk})^{-\frac{1}{2}}]\times exp[-({z(p_{j})}^{2}-2z(p_{j})z(p_{k})\delta_{jk}+{z(p_{k})}^{2}]/2(1-\delta_{jk}^{2})$ with $j=1,2,...,J-1;k=2,3...,J.$ Note that, on the condition $\sigma_{jk}$ is within the correlation range mentioned in Demirtas and Hedeker, (2011), the solution can be obtained uniquely. After generating the normal outcomes $(Z_{j})$, binary variables $(Y_{j})$ can be constructed by setting $Y_{j}=1$ if $Z_{j}\leq z(p_{j})$ and [math] if otherwise. Equivalently, by setting $Y_{j}=1$ if $Z_{j}\geq z(1-p_{j})$, binary variables could be created because there will be no change in the correlation. Note also that the discretization proportion corresponds to the proportion of 1’s which we observe by dichotomization.

Our main interest in this study is to generate an fBm via correlated random walk as an analogous method of generating Brownian motion through identically independently distributed random walk. For this purpose, the discretization context explained above can be implemented to the increments of fBm, in order to create a random walk with an explicit correlation structure. Therefore, we derive the correlated random walk by setting $Y_{j}=1$ if $Z_{j}\leq z(p_{j})$ and $Y_{j}=-1$ if $Z_{j}>z(p_{j})$. This creates the increments of correlated random walk and by taking the cumulative sum of them we can construct a correlated path. Notice that, Eqn. (3) is assured when binary variables are observed by setting $Y_{j}=1$ if $Z_{j}\leq z(p_{j})$ and [math] otherwise. Fortunately, binary variables $Y$ can also be assigned the values of $+1$ and $-1$ by taking the linear transformation and it can be checked that the link between the tetrachoric correlation and the phi-coefficient given in Eqn. (3) remains the same. Now, for simplicity let us consider the discretization proportions $p_{j}$ and $p_{k}$ equal the same proportion $p,$ in which case the link given in Eqn. (3) can be reduced to,

[TABLE]

where $z(p)$ represents $p$th quantile of standard normal distribution and

[TABLE]

is the cumulative distribution function of standard bivariate normal with correlation coefficient $\delta_{jk}$ where $f(z(p),z(p),\delta_{jk})=[2\pi^{-1}(1-\delta_{jk})^{-\frac{1}{2}}]\times exp[-({z(p)}^{2}-2z(p)^{2}\delta_{jk}+{z(p)}^{2}]/2(1-\delta_{jk}^{2})$ with $j=1,2,...,J-1;k=2,3...,J.$

The increment process, fGn, $\{B^{H}(n+1)-B^{H}(n):n=0,1,2,...\},$ are jointly Gaussian variates with zero mean and covariance given in Eqn. (2). Then, the correlation between $B^{H}(n+m+1)-B^{H}(n+m)$ and $B^{H}(n+1)-B^{H}(n)$ equals,

[TABLE]

with integer valued $m$ and Hurst parameter $H\in(0,1)$. As a consequence of these properties of fBm the one step dichotomized variates of these increments would have the following correlation, phi-coefficient,

[TABLE]

for $n=0,1,2,...$ by replacing the tetrachoric correlation with the correlation of fGn. And, here in this study, in order to generate fBm via correlated random walk, we first propose to generate correlated binary variates, $\pm 1$’s where $P(Y_{n}=1)=p,$ $P(Y_{n}=-1)=1-p$ with $m-$step correlation

[TABLE]

Second, we propose to take the cumulative sum of these binary random variates to construct a correlated trajectory. In such a simulation study, we have to answer two main questions which are,

• •

Does this simulated trajectory converge to fBm?

• •

If it does, for which values of probability $p$ is the convergence satisfied?

3 The convergence of correlated random walk

In this Section, we prove the convergence of the proposed correlated random walk given in Section 2 to fBm. Enriquez, (2004) introduces the construction of fBm using random walk with certain persistence parameters. The persistent random walk is defined as a discrete time process including the jumps of size $\pm 1,$ whose probability of making the same jump as the previous one is the parameter of persistence. Let $X^{\rho}$ be a discrete process with persistence parameter, $\rho\in[0,1]$ then,

1. i.

$X^{\rho}_{0}=0,$ $P(X^{\rho}_{1}=-1)=1/2$ and $P(X^{\rho}_{1}=+1)=1/2,$ 2. ii.

for all $n\geq 1,$ $\epsilon^{\rho}_{n}=X^{\rho}_{n}-X^{\rho}_{n-1}$ is equal to $-1$ or $+1$ a.s., 3. iii.

$P(\epsilon^{\rho}_{n+1}=\epsilon^{\rho}_{n}|\sigma(X_{k}^{\rho},0\leq k\leq\ n))=\rho$ for all $n\geq 1.$

Furthermore, by the help of conditioning on $\epsilon^{\rho}_{1},$ it can be easily seen that for all $n\geq 1,$ we have $P(\epsilon^{\rho}_{n}=\pm 1)=1/2.$ Hence, the correlation among $n$ distance time steps is,

[TABLE]

for all $n\geq 0,$ $m\geq 1$. In Enriquez, (2004) study additional randomness into persistence parameter $\rho$ is introduced where $P^{\rho}$ stand for the law of $X^{\rho}$ for a given random persistence. After forming a new probability measure $\mu$ on $[0,1]$, the law of persistent random walk corresponding to this $\mu$ is $P^{\mu}:\int_{0}^{1}{P^{\rho}d\mu(\rho)}$ and under $P^{\mu}$ for the process $X^{\mu}$, the correlation among $n$ distance time steps is computed by

[TABLE]

for all $n\geq 0$ , $m\geq 1$. Additionally, Enriquez, (2004) proves that the persistent random walk $X^{\mu}$ with the law of $P^{\mu}$ weakly converges to fBm $B_{H}(t)$ by Lemma 5.1 of Taqqu, (1975). The convergence can be shown in two steps. First, the summation of great number of trajectories converge to discrete Gaussian process by the Central Limit Theorem, and second this discrete Gaussian process converges to fBm when rescaled as it satisfies the correlation condition stated in Taqqu, (1975). Now, let us start by exploring the connection between the persistant random walk and the correlated binary variates we have proposed in Section 2.

Lemma 1

Let $\rho$ be the persistence parameter. Then,

[TABLE]

where $\delta^{H}_{1}=\frac{1}{2}[2^{2H}-2].$

**Proof: **Let $Y_{n}$ denote the jump size of random walk at the $n^{th}$ step. Then $\rho=P(Y^{\rho}_{n+1}=Y^{\rho}_{n}|\sigma(Y_{k}^{\rho},0\leq k\leq n))$ for all $n\geq 1.$ Therefore in terms of dichotomization we have

[TABLE]

for some $Z_{n}$ and $Z_{n+1}$ bivariate Normal random variables with correlation coefficient $\delta_{n,n+1}.$ Recall that, $P[Z_{n}\leq z(p)]=p$, $P[Z_{n}>z(p)]=1-p$ and

[TABLE]

Besides,

[TABLE]

and we obtain

[TABLE]

Hence, by Eqns.(10) and (12) we have

[TABLE]

Note that, by the arguments provided for obtaining Eqn. (6), we can replace $\delta_{n,n+1}$ with $\delta^{H}_{1}$ to conclude the statement. $\Box$

Proposition 1

Let $Y_{n}$ and $Y_{n+1}$ be the dichotomized variables of one step increments of fBm then the $n-$step correlation among them is

[TABLE]

**Proof: **Given $\sigma(Y_{k}^{\rho},0\leq k\leq n)$, the random variables $Y_{n}$ and $Y_{n+1}$ are the dichotomized variables of one step increments of fBm which are bivariate Gaussian. Therefore their correlation equals,

[TABLE]

and hence $E[Y_{n}Y_{n+1}]=4\Phi[z(p),z(p),\delta^{H}_{1}]-4p+1.$ On the other hand for all $n\geq 1$ we have

[TABLE]

and thus $E[Y_{n}Y_{n+1}]=2\rho-1=4\Phi[z(p),z(p),\delta^{H}_{1}]-4p+1.$ As a result, we have $\rho=2\Phi[z(p),z(p),\delta^{H}_{1}]-2p+1$ which coincides with the result obtained in Lemma 1. Consequently, as given in Proposition 1 of Enriquez, (2004), one can show that conditioning on $\sigma(Y_{k}^{\rho},0\leq k\leq m+n)$ by induction for all $m\geq 1,~{}n\geq 1,$ $E(Y_{m},Y_{m+n})=(2\rho-1)^{n}.$ Thus, by Lemma 1, in our proposed method, we take the correlation equal to

[TABLE]

$\Box$

As a result of these relations, it is observed that both the persistence parameter and the correlation coefficients can be written in terms of the parameter $p$ of the marginal distributions of the binary variates, $\pm 1$’s. In fact, note that the persistence parameter is a function of $p$ and $H,$ only. Due to the properties of fBm, the Hurst parameter H is a constant in the interval $[0,1],$ which therefore introducing randomness to the persistence parameter results in introducing randomness to the parameter $p.$

In Theorem 1, as a consequence of obtaining the link between the persistence parameter and the parameter of marginal proportions of the correlated binary variates, we are able to answer both questions introduced at the end of Section 2. We prove that the proposed correlated random walk with proportion satisfying an explicit density in fact converges to fBm with $H\in[1/2,1]$.

Theorem 1

For $H\in[\frac{1}{2},1]$ let $p$, the marginal parameter of correlated binary random variables $\pm 1,$ have the density

[TABLE]

for values of $p$ in the range $0\leq\frac{\left(4\Phi\left[z(p),z(p),\delta_{1}^{H}\right]-4p+1\right)^{n}-(2p-1)^{2}}{4p(1-p)}\leq 1$ and $(X^{H,i})_{i\geq 1}$ be a sequence of standardized independent process of this law then,

[TABLE]

where $a_{H}=\sqrt{\frac{H(2H-1)}{\Gamma{(3-2H)}}},$ N is the number of time steps and M be the number of trajectories, $B^{H}(t)$ is the fBm with Hurst index $H\in[1/2,1].$

**Proof: **In order for a persistent random walk to converge to fBm, the density imposed on the persistence parameter$\rho\in[\frac{1}{2},1]$ is given as;

[TABLE]

for values of $H\in[\frac{1}{2},1],$ see Enriquez, (2004). Now, let $Y_{k}^{H,i}$ be $k$ th step of $i$ th trajectory consisting of the correlated binary variables ($\pm 1$) with marginal proportions $p$ and $1-p$ respectively for $k\geq 1$ and $1\leq i\leq M$ where $M$ represents the number of trajectories. Now, let $X_{k}^{H,i}=X_{k-1}^{H,i}+Y_{k}^{H,i},k\geq 1$ with $X_{0}=0.$ Note that, $E(Y_{k}^{H,i})=2p-1,$ for $k\geq 1$ and so $E(X_{k}^{H,i})=k(2p-1).$ Also for all $k\geq 1$,

[TABLE]

By induction we have $E((X^{H,i}_{k})^{2})=2(2p-1)^{2}\sum_{j=2}^{k}(j-1)+k=k+k(k-1)(2p-1)^{2}$ for all $k\geq 1.$ Hence, $V(X^{H,i}_{k})=k(1-(2p-1)^{2})$ for all $k\geq 1.$

Then, by the central limit theorem

[TABLE]

is a discrete centered Gaussian process for $k\geq 1$ with $E[Z^{H}_{k}]=0$ and $V[Z^{H}_{k}]=\frac{V[X_{k}^{H}]}{4p(1-p)}$ since it is the sum of $M$ independent identically distributed random variables. Now for $k\geq 1$, let us define the increment process,

[TABLE]

which is the sum of independent, identically distributed random variables. Therefore, by the central limit theorem, it is Gaussian and stationary with $E[G^{H}_{k}]=0$ and $E[(G^{H}_{k})^{2}]=1$ for $k\geq 1.$ The n step correlation is;

[TABLE]

Now, the $n$ step correlation of $\pm 1^{\prime}$s, $\sigma_{k,k+n}$ is

[TABLE]

and a distribution can be imposed to $p$ by the help of the distribution of $\rho$ given in Eqn (19) for which the convergence is satisfied. Thus using the transformation of random variables on the density of $\rho,$ the density of $p$ is obtained as

[TABLE]

As a consequence, we observe that

[TABLE]

where $v=\frac{\left(\frac{\left(4\Phi\left[z(p),z(p),\delta_{1}^{H}\right]-4p+1\right)^{n}-(2p-1)^{2}}{4p(1-p)}\right)^{1/n}+1}{2}$ and $p_{l}$ and $p_{u}$ are the lower and upper values of $p$ which satisfy the inequality

[TABLE]

and $a^{2}_{H}$ satisfies

[TABLE]

as given in the proof of Thm 1 of Enriquez, (2004). Hence, by the Taqqu, (1975), Lemma 5.1 we conclude that the sum of the Bernoulli($p$) random variables stated in this Theorem converge to fBm with $H\in[1/2,1]$ $\Box$

4 A simple algorithm for generating fractional Brownian Motion using correlated random walk

The aim of this article is to propose an alternative yet simple algorithm which uses multivariate Bernoulli random variables to generate fBm through correlated random walk. In this paper, we mainly construct a random walk with the $n-$step correlation given in Proposition 1 and with a marginal parameter $p,$ satisfying an explicit density given in Theorem 1. We have proved that this correlated random walk converges weakly to fractional Brownian motion using the relation between $p$ and the persistence paramater $\rho$. This newly proposed algorithm generalizes the algorithm given in Enriquez, (2004) by carrying the information in the persistence parameter to the marginal parameter $p$ of the correlated Bernoulli random variables. It generalizes the construction given in Enriquez, (2004) since it uses the same persistence values but with the Bernoulli $p$ random variables as increments where $p$ has a certain density function, rather than using only Bernoulli $p=1/2$ random variables. The algorithm is proposed in a simple way to be followed and implemented by various groups of researchers working in different fields of applied science. This newly proposed algorithm has also theoretical interest since it is the exact analogous method of independent random walk converging to Brownian motion. In addition, in this algorithm one makes use of the correlation of the one step increments of fBm while generating correlated Bernoulli random variables.

Initially, we obtain the values of proportion $p$ which satisfy the density given in Eqn.(18). Notice that,

[TABLE]

where $v=\frac{\left(\frac{\left(4\Phi\left[z(p),z(p),\delta_{1}^{H}\right]-4p+1\right)^{n}-(2p-1)^{2}}{4p(1-p)}\right)^{1/n}+1}{2}.$ Equivalently,

[TABLE]

for $n=1.$ Therefore, the marginal parameter of the correlated Bernoulli random variables $p$ takes on the values which satisfy the Equation (25).

In the algorithm, by generating Uniform(0,1) random variable, $u,$ we first obtain the parameter $p$ of correlated Bernoulli random variables which ensure the equality (25). This implies that

[TABLE]

as given in Lemma 1 and Proposition 1. Now, consider $\xi(i-1)$ which is taken to be Bernoulli(p) random variable then note that

[TABLE]

have Bernoulli(p) distribution for all $i=1,...,n$ and $P(\xi(i)=\xi(i-1)|\sigma(\xi_{i},1\leq i\leq n-1))$ equals $2\Phi\left[z(p),z(p),\delta_{1}^{H}\right]-2p+1$ Hence, the $n$ step correlations among them satisfy $\frac{\left(4\Phi\left[z(p),z(p),\delta_{1}^{H}\right]-4p+1\right)^{n}-(2p-1)^{2}}{4p(1-p)}$ as requested. Converting these Bernoulli random variables to $\pm 1$s would not change the correlation coefficient among them. Thus, after having the Bernoulli(p) $\pm 1$ random variables with the desired correlation coefficients, the levels of the correlated walk is acquired by taking the summation of these increments. Finally, we use Theorem 1 which suggests that if the sums of these correlated random walk at each step which are scaled by the terms $M(Nt(2p-1))$ and $N^{H}\sqrt{M4p(1-p)}$ and multiplied with the constant $a_{H}$, the resulting process converges weakly to fBm as the number of trajectories $M$, and the number of steps $N$ become larger. Hence, to generate a trajectory of fBm, the steps to follow in this newly proposed algorithm are;

• •

Generate the value of $p$ from the distribution given in Theorem 1, that is find a value of $p$ which satisfy the equality (25).

• •

Generate

[TABLE]

and

[TABLE]

• •

Observe $\pm 1$ using

[TABLE]

• •

Repeat the first two steps $M$ times and observe $MXN$ multivariate Bernoulli ($\pm$1) random variables where the columns are generated with the marginal parameter $p,$ and the $n-$step correlation among the columns equals Proposition 1, (14).

• •

First take the cumulative sum of each row and form new columns using these cumulative sums. And then take the sum of resulting columns, subtract $M(Nt)(2p-1)$ from each column sum and divide by $N^{H}\sqrt{M4p(1-p)},$ and multiply the result by $a_{H}=\sqrt{\frac{H(2H-1)}{\Gamma{(3-2H)}}}$, for large values of $M$ and $N.$

The advantage of this algorithm is that, it is based on the construction which is done analogously to the construction of the standard Brownian motion generated by independent random walk. Therefore, it provides theoretical attraction. In addition, the correlated random walk proposed here still carries the information of one step correlation of the increments of fBm. The algorithm is easier to follow for all the researchers with only the knowledge of first year Statistics courses, which makes it easy to implement. It transforms the information of the persistence parameter to the marginal parameter of correlated Bernoulli random variables and instead of generating only Bernoulli($1/2$) with persistence parameter, in this algorithm we generate Bernoulli(p) random variables with a certain distribution assigned on $p.$ Hence this method generalizes the algorithm proposed in Enriquez, (2004), in the way that it uses the same persistence of increments of fBm but with Bernoulli(p) random variables with certain distribution assigned on $p$ rather than only using the Bernoulli($1/2$) variables. The number of computations and the memory required by this algorithm is same as Enriquez, (2004) algorithm. That is, for an $N$ step trajectory the number of computations is a power of $N$ between $1$ and $2$ and the memory required by the computations is a power of $N$ between [math] and $1,$ which makes the algorithm efficient compared to the other suggested methods. Comparisons can be found in Enriquez, (2004) study. The proposed algorithm is illustrated by Figure 1 with Hurst parameters, 0.55, 0.7 and 0.85 respectively. For each realization the Hurst parameters are estimated using Veitch and Abry, (1999) method and the results are presented.

In fact, the convergence is satisfied for less number of steps however the Veitch and Abry, (1999) method requires at least 10000 steps. For this reason, we have taken number of steps as 100000 in the above illustrations. We are also presenting a realization and its estimator for 10000 steps and 10000 paths in the following Figure 2.

We have also estimated the Hurst parameter of 30 realizations of our proposed algorithm using Veitch and Abry, (1999) Hurst index estimation method. The histogram of these estimators are illustrated in Figure 3 for realizations with H=0.7, 10000 steps, 10000 paths.