spatialstats.sample_empirical_variogram issue with scipy >= 1.13

[tayden]

The latest scipy version removes the _update() method from scipy.sparse.dok_matrix. This causes xdem.spatialstats.sample_empirical_variogram to throw an error, since it indirectly uses this method via skgstats.MetricSpace.

The version of scipy used by this package may need to be limited to <1.13 to avoid this.

xdem version: 0.0.18
OS: Ubuntu 22.04

Stack Trace

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[8], line 91
     85 plt.legend(loc="lower right")
     86 plt.show()
---> 91 df_vgm = xdem.spatialstats.sample_empirical_variogram(
     92     values=z_dh.data.squeeze(),
     93     gsd=dh.res[0],
     94     subsample=300,
     95     n_variograms=10,
     96     estimator="dowd",
     97     random_state=42,
     98 )
    100 func_sum_vgm, params_vgm = xdem.spatialstats.fit_sum_model_variogram(
    101     ["Gaussian", "Spherical"], empirical_variogram=df_vgm
    102 )
    103 xdem.spatialstats.plot_variogram(
    104     df_vgm,
    105     xscale_range_split=[100, 1000, 10000],
    106     list_fit_fun=[func_sum_vgm],
    107     list_fit_fun_label=["Standardized double-range variogram"],
    108 )

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/xdem/spatialstats.py:1498](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/xdem/spatialstats.py#line=1497), in sample_empirical_variogram(values, gsd, coords, subsample, subsample_method, n_variograms, n_jobs, verbose, random_state, **kwargs)
   1489     for i in range(n_variograms):
   1491         argdict = {
   1492             "i": i,
   1493             "imax": n_variograms,
   (...)
   1496             **kwargs,  # type: ignore
   1497         }
-> 1498         df_run = _wrapper_get_empirical_variogram(argdict=argdict)
   1500         list_df_run.append(df_run)
   1501 else:

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/xdem/spatialstats.py:1271](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/xdem/spatialstats.py#line=1270), in _wrapper_get_empirical_variogram(argdict)
   1267 argdict.pop("imax")
   1269 if argdict["subsample_method"] in ["cdist_equidistant", "cdist_point"]:
   1270     # Simple wrapper for the skgstat Variogram function for cdist methods
-> 1271     return _get_cdist_empirical_variogram(**argdict)
   1272 else:
   1273     # Aggregating several skgstat Variogram after iterative subsampling of specific points in the Raster
   1274     return _aggregate_pdist_empirical_variogram(**argdict)

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/xdem/spatialstats.py:1243](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/xdem/spatialstats.py#line=1242), in _get_cdist_empirical_variogram(values, coords, subsample_method, **kwargs)
   1241 # Filter corresponding arguments before passing to Variogram function
   1242 filtered_var_kwargs = {k: kwargs[k] for k in variogram_args if k in kwargs}
-> 1243 V = skg.Variogram(M, values=values, normalize=False, fit_method=None, **filtered_var_kwargs)
   1245 # Get bins, empirical variogram values, and bin count
   1246 bins, exp = V.get_empirical(bin_center=False)

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py:394](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py#line=393), in Variogram.__init__(self, coordinates, values, estimator, model, dist_func, bin_func, normalize, fit_method, fit_sigma, use_nugget, maxlag, samples, n_lags, verbose, **kwargs)
    392 fit_bounds = self._kwargs.get('fit_bounds') # returns None if empty
    393 fit_p0 = self._kwargs.get('fit_p0')
--> 394 self.fit(force=True, bounds=fit_bounds, p0=fit_p0)
    396 # finally check if any of the uncertainty propagation kwargs are set
    397 self._experimental_conf_interval = None

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py:1641](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py#line=1640), in Variogram.fit(self, force, method, sigma, bounds, p0, **kwargs)
   1638     self.cov = None
   1640 # if force, force a clean preprocessing
-> 1641 self.preprocessing(force=force)
   1643 # load the data
   1644 x = np.array(self.bins)

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py:1550](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py#line=1549), in Variogram.preprocessing(self, force)
   1529 """Preprocessing function
   1530 
   1531 Prepares all input data for the fit and transform functions. Namely,
   (...)
   1547 
   1548 """
   1549 # call the _calc functions
-> 1550 self._calc_diff(force=force)
   1551 self._calc_groups(force=force)

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py:1933](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py#line=1932), in Variogram._calc_diff(self, force)
   1930     return
   1932 # format into column-stack for faster calculation
-> 1933 diffs = self._format_values_stack(self.values)
   1935 # check if this is a cross-variogram
   1936 if self.is_cross_variogram:

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py:1879](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py#line=1878), in Variogram._format_values_stack(self, values)
   1871 """
   1872 Create a numpy column stack to calculate differences between two value arrays.
   1873 The format function will handle sparse matrices, as these do not include
   (...)
   1876 
   1877 """
   1878 # handle sparse matrix
-> 1879 if isinstance(self.distance_matrix, sparse.spmatrix):
   1880     # get triangular distance matrices
   1881     c = r = self.triangular_distance_matrix
   1883     # get the sparse CSR matrix

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py:1222](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/Variogram.py#line=1221), in Variogram.distance_matrix(self)
   1220 @property
   1221 def distance_matrix(self):
-> 1222     return self._X.dists

File [~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/MetricSpace.py:786](http://10.8.1.28:8888/lab/tree/Derek_big_bar/~/miniforge3/envs/derek/lib/python3.11/site-packages/skgstat/MetricSpace.py#line=785), in RasterEquidistantMetricSpace.dists(self)
    776 # remove possible duplicates (that would be summed by default)
    777 # from https://stackoverflow.com/questions/28677162/ignoring-duplicate-entries-in-sparse-matrix
    778 
   (...)
    783 # Solution 5+ times faster than the preceding, but relies on _update() which might change in scipy (which
    784 # only has an implemented method for summing duplicates, and not ignoring them yet)
    785 dok = sparse.dok_matrix((len(self.coords), len(self.coords)))
--> 786 dok._update(zip(zip(c, eq), d))
    787 dists = dok.tocsr()
    789 self._dists = dists

AttributeError: 'dok_matrix' object has no attribute '_update'

[tayden]

Scipy appears to be limited already on the main branch, but hasn't made it into a release yet.

Closing this, but I've submitted it to possibly help others find the solution until this fix is released.

[rhugonnet]

Thanks @tayden! Linked to #503 and mmaelicke/scikit-gstat#178

I will make a new release by the end of the week 😉, a couple other things need to be crammed in before