Estimate the Hurst exponent

The Hurst exponent $H \in (0, 1)$ measures the long memory of a series. Fractional Brownian motion $B^H$ has increments that are positively correlated for $H > \tfrac12$, negatively for $H < \tfrac12$, and independent at $H = \tfrac12$, where it is Brownian motion. This tutorial simulates fractional Brownian motion with a known $H$, reads $H$ back with four of the library's estimators, and says which estimator to trust for which data.

What you'll build

One path with $H = 0.7$, four estimates of $H$ from it — two of them with a standard error — and a rule for choosing among the eight estimators the library ships.

Prerequisites

• Rust: stochastic-rs, no features. Python: pip install stochastic-rs.
• The estimators live in stochastic_rs::stats::hurst (Higuchi and Variogram in stats::fractal_dim) and share one trait, HurstEstimator, whose estimate(path) returns a HurstResult with the point estimate, an optional standard error and the regression diagnostics.

Step 1 — simulate a path with a known H

Fractional Brownian motion is the Gaussian process with

Fbm builds the path from fractional Gaussian noise (Fgn), which it samples exactly by circulant embedding (Davies–Harte), so the path has this covariance up to floating-point rounding, not up to a discretisation error.

import stochastic_rs as srs

path = srs.PyFbm(0.7, 4096, t=1.0, seed=1).sample()  # 4096 points on [0, 1]

In Rust the same line is Fbm::<f64, _>::new(0.7, 4_096, Some(1.0), Deterministic::new(1)).sample(), which opens the example in step 2.

Step 2 — estimate H

// docs: tutorials/hurst-exponent
//! Estimates the Hurst exponent of a simulated fBm path with four estimators behind one trait.

use stochastic_rs::simd_rng::Deterministic;
use stochastic_rs::stats::fractal_dim::Higuchi;
use stochastic_rs::stats::hurst::Dfa;
use stochastic_rs::stats::hurst::Gph;
use stochastic_rs::stats::hurst::Wavelet;
use stochastic_rs::stochastic::process::fbm::Fbm;
use stochastic_rs::traits::HurstEstimator;
use stochastic_rs::traits::ProcessExt;

#[test]
fn estimate_hurst_of_an_fbm_path() {
  let path = Fbm::<f64, _>::new(0.7, 4_096, Some(1.0), Deterministic::new(1)).sample();

  // Every default here takes the level path, not its increments.
  let dfa = Dfa::default().estimate(path.view()).unwrap();
  let gph = Gph::default().estimate(path.view()).unwrap();
  let wavelet = Wavelet::default().estimate(path.view()).unwrap();
  let higuchi = Higuchi::default().estimate(path.view()).unwrap(); // H = 2 - D

  for h in [dfa.hurst, gph.hurst, wavelet.hurst, higuchi.hurst] {
    assert!((h - 0.7).abs() < 0.1, "H = {h}");
  }
  // Only Gph and Wavelet report an asymptotic standard error.
  assert!(gph.std_err.is_some() && wavelet.std_err.is_some() && dfa.std_err.is_none());
}

import stochastic_rs as srs

path = srs.PyFbm(0.7, 4096, t=1.0, seed=1).sample()

for est in (srs.Dfa(), srs.Gph(), srs.Wavelet()):
    r = est.estimate(path)
    se = "none" if r.std_err is None else f"{r.std_err:.3f}"
    print(f"{type(est).__name__:8s} H = {r.hurst:.3f}  std err = {se}")
print(f"Higuchi  H = {srs.Higuchi().estimate_hurst(path).hurst:.3f}")

# Dfa      H = 0.684  std err = none
# Gph      H = 0.634  std err = 0.053
# Wavelet  H = 0.720  std err = 0.032
# Higuchi  H = 0.690

Every default here takes the level path — the fBm itself — and differences or integrates it internally. Two things to know before feeding it your own data:

• Increments instead of levels. For a noise series such as fractional Gaussian noise or returns, pass take_differences=False to RescaledRange, Gph and Wavelet, and assume_integrated=True to Dfa only if you have already cumulated the noise into a walk.
• Higuchi returns a fractal dimension. Higuchi and Variogram estimate the fractal dimension $D$ of the curve; estimate_hurst (Python) or the HurstEstimator impl (Rust) converts it as $H = 2 - D$, clamped to $(0, 1)$. Call estimate for $D$ itself when a degenerate fit should stay visible instead of being squashed to a boundary value.

Step 3 — choose an estimator

Seven of the eight estimators measure the same thing, the self-similarity exponent of the series you hand them, by different routes; one measures something else.

Whittle (and FukasawaHurst in Python) implements the adapted Whittle estimator of Fukasawa, Takabatake and Westphal (2019): it estimates the Hurst exponent of a latent volatility process from log realised variance aggregated from intraday returns (FukasawaHurst.from_prices(closes), FukasawaHurst.from_log_rv(log_rv)). Handed a directly sampled path it does not fail — it returns a number unrelated to that path's exponent — so build its input from prices, not from a simulated Fgn.

The estimators also have hard minimum lengths with their defaults: 32 observations for RescaledRange and Dfa, 64 for Gph and Wavelet, 30 for Whittle. Those are validation floors, not sample sizes to aim for — a log-log regression wants thousands of points.

Result

On this path the four estimates land within 0.07 of the true 0.7, and the two standard errors (0.053 for GPH, 0.032 for the wavelet estimator) say how far apart to expect estimates on paths of this length. Rerun with other seeds and lengths and the spread of each estimator is the thing to watch: that spread, more than any one estimate, is what the choice of estimator changes.

Where to go next

• Choosing a Hurst estimator on the statistics page, with the per-estimator notes and their references
• Fractional Brownian motion and the other fractional processes, which all take hurst as a parameter
• GPU paths on a free Colab GPU — fractional Gaussian noise has its own device pipeline

References

• Hurst, H. E. (1951). Long-term storage capacity of reservoirs. doi:10.1061/TACEAT.0006518
• Geweke, J. and Porter-Hudak, S. (1983). The estimation and application of long memory time series models. doi:10.1111/j.1467-9892.1983.tb00371.x
• Higuchi, T. (1988). Approach to an irregular time series on the basis of the fractal theory. doi:10.1016/0167-2789(88)90081-4
• Peng, C.-K., Buldyrev, S. V., Havlin, S., Simons, M., Stanley, H. E. and Goldberger, A. L. (1994). Mosaic organization of DNA nucleotides. doi:10.1103/PhysRevE.49.1685
• Veitch, D. and Abry, P. (1999). A wavelet-based joint estimator of the parameters of long-range dependence. doi:10.1109/18.761330
• Fukasawa, M., Takabatake, T. and Westphal, R. (2019). Is volatility rough?. arXiv:1905.04852