.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/02_cov_model/07_roughness.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_02_cov_model_07_roughness.py: Roughness ========= The roughness describes the power-law behavior of a variogram at the origin ([Wu2016]_): .. math:: \gamma(r) \sim c \cdot r^\alpha \quad (r \to 0) The exponent :math:`\alpha` is the roughness information, bounded by :math:`0 \le \alpha \le 2`. A value of 0 corresponds to a nugget effect and 2 is the smooth limit for random fields. You can access it via :any:`CovModel.roughness`. On a log-log plot, the slope of the variogram near the origin approaches this value. .. GENERATED FROM PYTHON SOURCE LINES 17-23 .. code-block:: Python import numpy as np import matplotlib.pyplot as plt import gstools as gs .. GENERATED FROM PYTHON SOURCE LINES 24-28 Variogram behavior near the origin ---------------------------------- Compare variograms near the origin for models with different roughness. .. GENERATED FROM PYTHON SOURCE LINES 28-35 .. code-block:: Python models = [ gs.Exponential(len_scale=1.0), gs.Gaussian(len_scale=1.0), gs.Stable(len_scale=1.0, alpha=0.7), ] .. GENERATED FROM PYTHON SOURCE LINES 36-38 Use a small-lag grid and fit the slope on a log-log scale to estimate the roughness numerically. .. GENERATED FROM PYTHON SOURCE LINES 38-55 .. code-block:: Python r = np.logspace(-4, -1, 200) fig, ax = plt.subplots() for model in models: gamma = model.variogram(r) slope = np.polyfit(np.log(r[:20]), np.log(gamma[:20]), 1)[0] ax.loglog( r, gamma, label=rf"{model.name} ($\alpha={model.roughness:.2f}$)" ) print(f"{model.name}: roughness={model.roughness:.2f}, fit={slope:.2f}") ax.set_xlabel("r") ax.set_ylabel(r"$\gamma(r)$") ax.set_title("Variogram near the origin (log-log)") ax.legend() .. image-sg:: /examples/02_cov_model/images/sphx_glr_07_roughness_001.png :alt: Variogram near the origin (log-log) :srcset: /examples/02_cov_model/images/sphx_glr_07_roughness_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Exponential: roughness=1.00, fit=1.00 Gaussian: roughness=2.00, fit=2.00 Stable: roughness=0.70, fit=0.70 .. GENERATED FROM PYTHON SOURCE LINES 56-57 A nugget masks the power-law behavior at the origin, so roughness is 0. .. GENERATED FROM PYTHON SOURCE LINES 57-61 .. code-block:: Python nugget_model = gs.Gaussian(nugget=1.0) print("Gaussian with nugget roughness:", nugget_model.roughness) .. rst-class:: sphx-glr-script-out .. code-block:: none Gaussian with nugget roughness: 0.0 .. GENERATED FROM PYTHON SOURCE LINES 62-76 Roughness and random fields --------------------------- From the theory of fractal stochastic processes (Mandelbrot and Van Ness 1968) ([Mandelbrot1968]_), the roughness can be interpreted as: 1. Persistent (:math:`1 < \alpha \le 2`): smooth behavior, slowly increasing variogram, long memory (e.g. Gaussian-like). 2. Antipersistent (:math:`0 < \alpha < 1`): rough behavior, steep variogram near the origin, short memory. 3. No memory (:math:`\alpha = 1`): linear slope at the origin (Exponential). The Integral model lets us control the roughness with its parameter ``nu``. For this model, the roughness is given by :math:`\alpha = \min(2, \nu)`. .. GENERATED FROM PYTHON SOURCE LINES 78-79 Set up a common grid and plotting scales so the realizations are comparable. .. GENERATED FROM PYTHON SOURCE LINES 79-87 .. code-block:: Python sep = 100 # separation point (mid field) ext = 10 # field extend imext = 2 * [0, ext] x = y = np.linspace(0, ext, 2 * sep + 1) xm = np.linspace(0, 5, 1001) vmin, vmax = -3, 3 .. GENERATED FROM PYTHON SOURCE LINES 88-89 Create Integral models with the same integral scale but different roughness. .. GENERATED FROM PYTHON SOURCE LINES 89-93 .. code-block:: Python rough = [0.1, 1, 10] mod = [gs.Integral(dim=2, integral_scale=1, nu=nu) for nu in rough] .. GENERATED FROM PYTHON SOURCE LINES 94-98 .. note:: Rough fields require a higher ``mode_no`` so the randomization method samples the high-frequency part of the spectrum sufficiently well. .. GENERATED FROM PYTHON SOURCE LINES 100-101 Generate random field realizations with a shared seed for fair comparison. .. GENERATED FROM PYTHON SOURCE LINES 101-109 .. code-block:: Python srf = [] for m in mod: mode_no = int(1000 / m.roughness) # increase mode_no for rough fields s = gs.SRF(m, seed=20110101, mode_no=mode_no) s.structured((x, y)) srf.append(s) .. GENERATED FROM PYTHON SOURCE LINES 110-111 Plot variograms, fields, and cross sections column-wise by roughness. .. GENERATED FROM PYTHON SOURCE LINES 111-142 .. code-block:: Python fig, axes = plt.subplots(3, 3, figsize=(10, 10)) for i in range(3): label = rf"$\alpha(\gamma)={mod[i].roughness:.2f}$" axes[0, i].plot(xm, mod[i].variogram(xm), label=label) axes[0, i].legend(loc="lower right") axes[0, i].set_ylim([-0.05, 1.05]) axes[0, i].set_xlim([-0.25, 5.25]) axes[0, i].grid() axes[1, i].imshow( srf[i].field.T, origin="lower", interpolation="bicubic", vmin=vmin, vmax=vmax, extent=imext, ) axes[1, i].axhline(y=5, color="k") axes[2, i].plot(x, srf[i].field[:, sep]) axes[2, i].set_ylim([vmin, vmax]) axes[2, i].set_xlim([0, ext]) axes[2, i].grid() axes[0, 0].set_ylabel(r"$\gamma$") axes[1, 0].set_ylabel(r"$y$") axes[2, 0].set_ylabel(r"$f$") axes[2, 0].set_xlabel(r"$x$") axes[2, 1].set_xlabel(r"$x$") axes[2, 2].set_xlabel(r"$x$") fig.tight_layout() .. image-sg:: /examples/02_cov_model/images/sphx_glr_07_roughness_002.png :alt: 07 roughness :srcset: /examples/02_cov_model/images/sphx_glr_07_roughness_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 144-153 Illustration of the impact of the roughness information on spatial random fields. Each column shows the variogram on top, a single field realization in the middle and the highlighted cross section of the field on the bottom. The left column shows a very rough (antipersistent) field with a steep increase of the variogram at the origin, the middle column shows an Exponential like variogram with a linear slope at the origin resulting in a less rough (no memory) field and the right column shows a Gaussian like variogram together with a very smooth (persistent) field. All variograms have the same integral scale and x- and y-axis are given in multiples of the integral scale. .. GENERATED FROM PYTHON SOURCE LINES 155-163 References ---------- .. [Wu2016] Wu, W.-Y., and C. Y. Lim. 2016. "ESTIMATION OF SMOOTHNESS OF A STATIONARY GAUSSIAN RANDOM FIELD." Statistica Sinica 26 (4): 1729-1745. https://doi.org/10.5705/ss.202014.0109. .. [Mandelbrot1968] Mandelbrot, B. B., and J. W. Van Ness. 1968. "Fractional Brownian Motions, Fractional Noises and Applications." SIAM Review 10 (4): 422-437. https://doi.org/10.1137/1010093. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 17.181 seconds) .. _sphx_glr_download_examples_02_cov_model_07_roughness.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 07_roughness.ipynb <07_roughness.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 07_roughness.py <07_roughness.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 07_roughness.zip <07_roughness.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_