From ff51a5360818418a495277b90cac92ec1fe505e1 Mon Sep 17 00:00:00 2001 From: "selman.ozleyen" Date: Wed, 1 Oct 2025 13:34:15 +0200 Subject: [PATCH 1/4] add numerical stability --- src/squidpy/gr/_sepal.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/src/squidpy/gr/_sepal.py b/src/squidpy/gr/_sepal.py index 95d740997..fe386cff5 100644 --- a/src/squidpy/gr/_sepal.py +++ b/src/squidpy/gr/_sepal.py @@ -290,8 +290,11 @@ def _entropy( """Get entropy of an array.""" xnz = xx[xx > 0] xs: np.float64 = np.sum(xnz) + if xs == 0: + return 0.0 + eps = np.finfo(np.float64).eps # ~2.22e-16 xn = xnz / xs - xl = np.log(xn) + xl = np.log(np.maximum(xn, eps)) return float((-xl * xn).sum()) From fbc31ad20a9f23f484742180fd0d0004a044cae7 Mon Sep 17 00:00:00 2001 From: "selman.ozleyen" Date: Wed, 1 Oct 2025 13:39:53 +0200 Subject: [PATCH 2/4] simplify unnecessary parameters and arrays --- src/squidpy/gr/_sepal.py | 17 ++++------------- 1 file changed, 4 insertions(+), 13 deletions(-) diff --git a/src/squidpy/gr/_sepal.py b/src/squidpy/gr/_sepal.py index fe386cff5..c645f964b 100644 --- a/src/squidpy/gr/_sepal.py +++ b/src/squidpy/gr/_sepal.py @@ -202,14 +202,13 @@ def _score_helper( @njit(fastmath=True) def _diffusion( conc: NDArrayA, - laplacian: Callable[[NDArrayA, NDArrayA, NDArrayA], float], + laplacian: Callable[[NDArrayA, NDArrayA], float], n_iter: int, sat: NDArrayA, sat_idx: NDArrayA, unsat: NDArrayA, unsat_idx: NDArrayA, dt: float = 0.001, - D: float = 1.0, thresh: float = 1e-8, ) -> float: """Simulate diffusion process on a regular graph.""" @@ -217,15 +216,14 @@ def _diffusion( entropy_arr = np.zeros(n_iter) prev_ent = 1.0 nhood = np.zeros(sat_shape) - weights = np.ones(sat_shape) for i in range(n_iter): for j in range(sat_shape): nhood[j] = np.sum(conc[sat_idx[j]]) - d2 = laplacian(conc[sat], nhood, weights) + d2 = laplacian(conc[sat], nhood) dcdt = np.zeros(conc_shape) - dcdt[sat] = D * d2 + dcdt[sat] = d2 conc[sat] += dcdt[sat] * dt conc[unsat] += dcdt[unsat_idx] * dt # set values below zero to 0 @@ -246,7 +244,6 @@ def _diffusion( def _laplacian_rect( centers: NDArrayA, nbrs: NDArrayA, - h: float, ) -> NDArrayA: """ Five point stencil approximation on rectilinear grid. @@ -254,8 +251,6 @@ def _laplacian_rect( See `Wikipedia `_ for more information. """ d2f: NDArrayA = nbrs - 4 * centers - d2f = d2f / h**2 - return d2f @@ -264,7 +259,6 @@ def _laplacian_rect( def _laplacian_hex( centers: NDArrayA, nbrs: NDArrayA, - h: float, ) -> NDArrayA: """ Seven point stencil approximation on hexagonal grid. @@ -275,10 +269,7 @@ def _laplacian_hex( Curtis D. Benster, L.V. Kantorovich, V.I. Krylov, ISBN-13: 978-0486821603. """ - d2f: NDArrayA = nbrs - 6 * centers - d2f = d2f / h**2 - d2f = (d2f * 2) / 3 - + d2f: NDArrayA = (2.0 * nbrs - 12.0 * centers) / 3.0 return d2f From 901cc0f608bb2901d03bbf4c150007aa577ddba7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Selman=20=C3=96zleyen?= <32667648+selmanozleyen@users.noreply.github.com> Date: Wed, 15 Oct 2025 19:06:29 +0200 Subject: [PATCH 3/4] Update src/squidpy/gr/_sepal.py Co-authored-by: Philipp A. --- src/squidpy/gr/_sepal.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/squidpy/gr/_sepal.py b/src/squidpy/gr/_sepal.py index c645f964b..edec28331 100644 --- a/src/squidpy/gr/_sepal.py +++ b/src/squidpy/gr/_sepal.py @@ -278,7 +278,7 @@ def _laplacian_hex( def _entropy( xx: NDArrayA, ) -> float: - """Get entropy of an array.""" + """Compute Shannon entropy of an array of probability values (in nats).""" xnz = xx[xx > 0] xs: np.float64 = np.sum(xnz) if xs == 0: From 3f2d4c2726c509ab0c6bfa3024e685eb2b16723b Mon Sep 17 00:00:00 2001 From: selmanozleyen Date: Wed, 15 Oct 2025 19:24:44 +0200 Subject: [PATCH 4/4] add documentation --- src/squidpy/gr/_sepal.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/src/squidpy/gr/_sepal.py b/src/squidpy/gr/_sepal.py index edec28331..a5260050c 100644 --- a/src/squidpy/gr/_sepal.py +++ b/src/squidpy/gr/_sepal.py @@ -281,9 +281,13 @@ def _entropy( """Compute Shannon entropy of an array of probability values (in nats).""" xnz = xx[xx > 0] xs: np.float64 = np.sum(xnz) - if xs == 0: - return 0.0 eps = np.finfo(np.float64).eps # ~2.22e-16 + if xs < eps: + # 0 because + # xn represents probabilities + # and p(x)=0 is taken as 0 entropy + # see https://stats.stackexchange.com/a/433096 + return 0.0 xn = xnz / xs xl = np.log(np.maximum(xn, eps)) return float((-xl * xn).sum())