Merge 2f3145e into 63c3da4

docs/source/coregistration.rst

@@ -153,7 +153,7 @@ ICP
 - **Does not support weights**
 - **Recommended for:** Data with low noise and a high relative rotation.
 
-Iterative Closest Point (ICP) coregistration works by iteratively moving the data until it fits the reference as well as possible.
+Iterative Closest Point (ICP) coregistration, which is based on `Besl and McKay (1992) <https://doi.org/10.1117/12.57955>`_, works by iteratively moving the data until it fits the reference as well as possible.
 The DEMs are read as point clouds; collections of points with X/Y/Z coordinates, and a nearest neighbour analysis is made between the reference and the data to be aligned.
 After the distances are calculated, a rigid transform is estimated to minimise them.
 The transform is attempted, and then distances are calculated again.

docs/source/intro_accuracy_precision.rst

@@ -22,17 +22,33 @@ calculations consistent, reproducible, and easy.
 Accuracy and precision
 ----------------------
 
-Both `accuracy and precision <https://en.wikipedia.org/wiki/Accuracy_and_precision>`_ are important factors to account
-for when analyzing DEMs:
+`Accuracy and precision <https://en.wikipedia.org/wiki/Accuracy_and_precision>`_ describe random and systematic errors,
+respectively.
 
-- the **accuracy** (systematic error) of a DEM describes how close a DEM is to the true location of measured elevations on the Earth's surface,
-- the **precision** (random error) of a DEM describes the typical spread of its error in measurement, independently of a possible bias from the true positioning.
+*Note: sometimes "accuracy" is also used to describe both types of errors, and "trueness" systematic errors, as defined
+in* `ISO 5725-1 <https://www.iso.org/obp/ui/#iso:std:iso:5725:-1:ed-1:v1:en>`_ *. Here, we used accuracy for systematic
+errors as, to our knowledge, it is a more commonly used terminology in remote sensing applications.*
 
 .. figure:: imgs/precision_accuracy.png
     :width: 80%
 
-Source: `antarcticglaciers.org <http://www.antarcticglaciers.org/glacial-geology/dating-glacial-sediments2/precision-
-and-accuracy-glacial-geology/>`_, accessed 29.06.21.
+    Source: `antarcticglaciers.org <http://www.antarcticglaciers.org/glacial-geology/dating-glacial-sediments2/precision-and-accuracy-glacial-geology/>`_, accessed 29.06.21.
+
+
+For DEMs, we thus have:
+
+- **DEM accuracy** (systematic error) describes how close a DEM is to the true location of measured elevations on the Earth's surface,
+- **DEM precision** (random error) of a DEM describes the typical spread of its error in measurement, independently of a possible bias from the true positioning.
+
+The spatial structure of DEMs complexifies the notion of accuracy and precision, however. Spatially structured
+systematic errors are often related to the gridded nature of DEMs, creating **affine biases** while other, **specific
+biases** exist at the pixel scale. For random errors, a variability in error magnitude or **heteroscedasticity** exists
+across the DEM, while spatially structured patterns of errors are linked to **spatial correlations**.
+
+.. figure:: https://github.com/rhugonnet/dem_error_study/blob/main/figures/fig_2.png?raw=true
+    :width: 100%
+
+    Source: `Hugonnet et al. (2022) <https://doi.org/10.1109/jstars.2022.3188922>`_.
 
 Absolute or relative accuracy
 -----------------------------
@@ -47,6 +63,8 @@ TODO: Add another little schematic!
 Optimizing DEM absolute accuracy
 --------------------------------
 
+
+
 Shifts due to poor absolute accuracy are common in elevation datasets, and can be easily corrected by performing a DEM
 co-registration to precise and accurate, quality-controlled elevation data such as `ICESat <https://icesat.gsfc.nasa.
 gov/icesat/>`_ and `ICESat-2 <https://icesat-2.gsfc.nasa.gov/>`_.
@@ -96,6 +114,6 @@ The tools for quantifying DEM precision are described in :ref:`spatialstats`.
     Functions that are used in several examples create duplicate examples intead of being merged into the list.
     Circumventing manually by selecting functions used only once in each example for now.
 
-.. minigallery:: xdem.spatialstats.neff_circ xdem.spatialstats.plot_1d_binning
+.. minigallery:: xdem.spatialstats.infer_heteroscedasticity_from_stable xdem.spatialstats.get_variogram_model_func xdem.spatialstats.sample_empirical_variogram
     :add-heading: Examples that use spatial statistics functions
 

docs/source/spatialstats.rst

@@ -6,17 +6,16 @@ Spatial statistics
 Spatial statistics, also referred to as `geostatistics <https://en.wikipedia.org/wiki/Geostatistics>`_, are essential
 for the analysis of observations distributed in space.
 To analyze DEMs, xDEM integrates spatial statistics tools specific to DEMs described in recent literature,
-in particular in `Rolstad et al. (2009) <https://doi.org/10.3189/002214309789470950>`_,
-`Dehecq et al. (2020) <https://doi.org/10.3389/feart.2020.566802>`_ and
-`Hugonnet et al. (2021) <https://doi.org/10.1038/s41586-021-03436-z>`_. The implementation of these methods relies
-partly on the package `scikit-gstat <https://mmaelicke.github.io/scikit-gstat/index.html>`_.
+in particular in `Hugonnet et al. (2022) <https://doi.org/10.1109/jstars.2022.3188922>`_ and
+`Rolstad et al. (2009) <https://doi.org/10.3189/002214309789470950>`_. The implementation of these methods relies
+partially on the package `scikit-gstat <https://mmaelicke.github.io/scikit-gstat/index.html>`_.
 
 The spatial statistics tools can be used to assess the precision of DEMs (see the definition of precision in :ref:`intro`).
 In particular, these tools help to:
 
-    - account for non-stationarities of elevation measurement errors (e.g., varying precision of DEMs with terrain slope),
-    - quantify the spatial correlation of measurement errors in DEMs (e.g., native spatial resolution, instrument noise),
-    - estimate robust errors for observations integrated in space (e.g., average or sum of samples),
+    - account for elevation heteroscedasticity (e.g., varying precision with terrain slope),
+    - quantify the spatial correlation of errors in DEMs (e.g., native spatial resolution, instrument noise),
+    - estimate robust errors for observations analyzed in space (e.g., average or sum of elevation, or of elevation changes),
     - propagate errors between spatial ensembles at different scales (e.g., sum of glacier volume changes).
 
 .. contents:: Contents 
@@ -73,7 +72,7 @@ Due to the sparsity of synchronous acquisitions, elevation data cannot be easily
 times. Thus, stable terrain is used a proxy to assess the precision of a DEM on all its terrain,
 including moving terrain that is generally of greater interest for analysis.
 
-As shown in Hugonnet et al. (in prep), accounting for :ref:`spatialstats_nonstationarity` is needed to reliably
+As shown in `Hugonnet et al. (2022) <https://doi.org/10.1109/jstars.2022.3188922>`_, accounting for :ref:`spatialstats_heterosc` is needed to reliably
 use stable terrain as a proxy for other types of terrain.
 
 .. _spatialstats_metrics:
@@ -99,7 +98,7 @@ To estimate the pixel-wise measurement error for elevation data, two issues aris
     1. The dispersion :math:`\sigma_{dh}` cannot be estimated directly on changing terrain,
     2. The dispersion :math:`\sigma_{dh}` can show important non-stationarities.
 
-The section :ref:`spatialstats_nonstationarity` describes how to quantify the measurement error as a function of
+The section :ref:`spatialstats_heterosc` describes how to quantify the measurement error as a function of
 several explanatory variables by using stable terrain as a proxy.
 
 Spatially-integrated elevation measurement error
@@ -129,23 +128,21 @@ and use those to integrate and propagate measurement errors in space.
 Workflow for DEM precision estimation
 -------------------------------------
 
-.. _spatialstats_nonstationarity:
+.. _spatialstats_heterosc:
 
-Non-stationarity in elevation measurement errors
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+Elevation heteroscedasticity
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
-Elevation data contains significant non-stationarities in elevation measurement errors.
+Elevation data contains significant variability in measurement errors.
 
-xDEM provides tools to **quantify** these non-stationarities along several explanatory variables,
-**model** those numerically to estimate an elevation measurement error, and **standardize** them for further analysis.
+xDEM provides tools to **quantify** this variability using explanatory variables, **model** those numerically to
+estimate a function predicting elevation error, and **standardize** data for further analysis.
 
-Quantify and model non-stationarites
-""""""""""""""""""""""""""""""""""""
+Quantify and model heteroscedasticity
+"""""""""""""""""""""""""""""""""""""
 
-Non-stationarities in elevation measurement errors correspond to a variability of the precision in the elevation
-observations with certain explanatory variables that can be terrain- or instrument-related.
-In statistical terms, it corresponds to an `heteroscedasticity <https://en.wikipedia.org/wiki/Heteroscedasticity>`_
-of elevation observations.
+Elevation `heteroscedasticity <https://en.wikipedia.org/wiki/Heteroscedasticity>`_ corresponds to a variability in
+precision of elevation observations, that are linked to terrain or instrument variables.
 
 .. math::
     \sigma_{dh} = \sigma_{dh}(\textrm{var}_{1},\textrm{var}_{2}, \textrm{...}) \neq \textrm{constant}
@@ -158,7 +155,7 @@ Owing to the large number of samples of elevation data, we can easily estimate t
         :lines: 18-19
         :language: python
 
-.. plot:: code/spatialstats_nonstationarity_slope.py
+.. plot:: code/spatialstats_heterosc_slope.py
     :width: 90%
 
 The most common explanatory variables are:
@@ -167,7 +164,7 @@ The most common explanatory variables are:
     - the quality of stereo-correlation that can explain a large part of the measurement error of DEMs generated by stereophotogrammetry,
     - the interferometric coherence that can explain a large part of the measurement error of DEMs generated by `InSAR <https://en.wikipedia.org/wiki/Interferometric_synthetic-aperture_radar>`_.
 
-Once quantified, the non-stationarities can be modelled numerically by linear interpolation across several
+Once quantified, elevation heteroscedasticity can be modelled numerically by linear interpolation across several
 variables using :func:`xdem.spatialstats.interp_nd_binning`.
 
 .. literalinclude:: code/spatialstats.py
@@ -219,8 +216,8 @@ average measurement error of the pixels in the subsample, evaluated through the
 Estimating the standard error of the mean of the standardized data :math:`\sigma_{\overline{z_{dh}}}\vert_{\mathbb{S}}`
 requires an analysis of spatial correlation and a spatial integration of this correlation, described in the next sections.
 
-.. minigallery:: xdem.spatialstats.nd_binning
-        :add-heading: Examples that deal with non-stationarities
+.. minigallery:: xdem.spatialstats.infer_heteroscedasticity_from_stable xdem.spatialstats.nd_binning
+        :add-heading: Examples that deal with elevation heteroscedasticity
         :heading-level: "
 
 .. _spatialstats_corr:
@@ -310,7 +307,7 @@ This can be performed through the function :func:`xdem.spatialstats.fit_sum_mode
         :lines: 31
         :language: python
 
-.. minigallery:: xdem.spatialstats.sample_empirical_variogram
+.. minigallery:: xdem.spatialstats.infer_spatial_correlation_from_stable xdem.spatialstats.sample_empirical_variogram
         :add-heading: Examples that deal with spatial correlations
         :heading-level: "