# Computing the Hurst exponent or fractal dimension of fractional Brownian motion

The Hurst exponent is related to the fractal dimension by noticing that the fractal dimension $D$ is equal to $2-H$, where $d$ is the intrinsic dimension and $H$ is the Hurst exponent, for 1-D fractional Brownian motion.

I can generate an fBm series in Mathematica as follows:

tlow = 1;
thigh = 1000;
tinc = 1;
hurst = 0.4;
dataz=RandomFunction[FractionalBrownianMotionProcess[hurst], {tlow, thigh, tinc}, 1];


However, using the method described here to compute the fractal dimension, I am unable to get close to recovering the Hurst exponent used to generate the fBm.

For example:

stdev[{x__}] := StandardDeviation[{x}]
stdev[{x_}] := 1
dataf = RotateRight[Mean[R /@ #/stdev /@ #] & /@ (Partition[dataz, 2^#] & /@ Range)]
Fit[Log[dataf], {1, x}, x]


Is this a numerical error or are there better ways to compute Hurst exponents and fractal dimension in Mathematica?

• What is R[] ? – Dr. belisarius Aug 31 '14 at 4:08
• I edited your code for readability. Please undo the edit if you disagree – Dr. belisarius Aug 31 '14 at 4:29
• I think he uses R as Rescaled Range (The rescaled range is a statistical measure of the variability of a time series...) Rescaled Range. Also see: Hust Exponent. I think I've got a notebook on this I'll try to dig it out of my archive. – Jagra Aug 31 '14 at 5:17
• You have some programs in Mathematica to compute fractal dimensions, Hurst coefficient and multi fractals for time series and images in the book Fractal Geography by Andre Dauphine (2013) – user19457 Aug 31 '14 at 7:14

Mathematica's estimation routines are able to recover the Hurst exponent from the sample:

BlockRandom[SeedRandom["mathematica.SE/58539"];
tlow = 1; thigh = 1000; tinc = 1; hurst = 0.4;
dataz = RandomFunction[FractionalBrownianMotionProcess[hurst], {tlow, thigh, tinc}, 1]];

FindProcessParameters[dataz, FractionalBrownianMotionProcess[h]]


{h -> 0.397063}

so the issue must be with your code.

This late answer gives some details of how to calculate a dimension for fractional Brownian motion. The documentation at Help > FractionalBrownianMotionProcess > Basic Examples shows that the variance

 Variance[FractionalBrownianMotionProcess[\[Mu], \[Sigma], h][t]]


is $t^{2 h} \sigma^2$, which is the same as that on page 84 of The Science of Fractal Images edited by Peitgen and Saupe. They write $\rm{Var}[X(t_2)-X(t_1)]=\sigma^2~Abs[t_2-t_1]^{2~h}$, where the Hurst exponent satisfies $0<h<1$, and the fractal dimension is $D=2-h$.

The strategy is to find the variance $V$ of data points separated by increasingly large time differences $k=t_2-t_1$. The slope and intercept of the best fit line on a log-log plot of $V$ versus $k$ give the Hurst exponent $h$ and volatility $\sigma$, respectively. Thus, $Log[V] = 2 \log[\sigma] + 2 h \log[k]$. The question assumes $\mu=0$ and $\sigma=1$.

Block[{tlow, thigh, tinc, hurst, dataz, u, k, fit},
SeedRandom["mathematica.SE/58539"];
tlow = 1;
thigh = 1000;
tinc = 1;
hurst = 0.4;
dataz = RandomFunction[
FractionalBrownianMotionProcess[hurst],
{tlow, thigh, tinc}, 1];
u = Table[
{N[Log[k]], Log[Variance[Differences[dataz, 1, k]]]},
{k, 1, 500}];

fit = Fit[u, {1, k}, k];

ListLogLogPlot[u,
Frame -> True,
FrameLabel -> {"Log[k] = Log[t2-t1]", "Log[Variance]"},
PlotLabel -> "Log[V] = " <> ToString[fit],
Epilog -> {
Text["\[Sigma] = " <> ToString[E^(fit[]/2)], {0.0, 1.3}],
Text["h = " <> ToString[fit[[2, 1]]/2], {0.0, 1.0}],
Thick, Red,
Line[{Log[{1, fit /. k -> 1}], Log[{6, fit /. k -> 6}]}]},
BaseStyle -> {FontSize -> 14}]
] 