docs: organise by Diátaxis mode, and keep the old URLs working

CHANGELOG.md

@@ -14,7 +14,7 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
   and
   the last training report. No pickle on either side, so loading one cannot execute code.
 
-  A reloaded map **continues the same random stream**, so training it further gives what never
+  A reloaded map continues the same random stream, so training it further gives what never
   stopping would have given. Both the seed and the generator's current state are saved: the seed
   alone
   would restart the stream and a resumed run would quietly diverge from an uninterrupted one.
@@ -31,13 +31,13 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
   and
   the member list, JSON, format version and weight shape are all validated before a map is built.
   What
-  that buys is "cannot execute code", not "safe to load anything" — an `.npz` is a zip, so a hostile
+  that buys is "cannot execute code", not "safe to load anything". An `.npz` is a zip, so a hostile
   file can still attempt resource exhaustion. `docs/save-and-load.md` states both halves.
 
 - **`SOM.last_report`**, a frozen `TrainingReport`: mode, iterations, sample count, seed, the final
   learning rate and radius, quantization error, the `python_som` and NumPy versions, and wall time.
   `train()` still returns the quantization error, so nothing that used the return value changes.
-  `final_learning_rate` is `None` for batch training, which has no step size — Eq. (8) is a weighted
+  `final_learning_rate` is `None` for batch training, which has no step size: Eq. (8) is a weighted
   mean, and reporting the unused initial value would read as though it had been applied.
 
 - **`SOMConfig`**, the serialisable description of a map, available as `som.config()`.
@@ -85,7 +85,7 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 - The package is now a pure functional core with a thin shell around it. `python_som._core` holds
   every numeric decision as functions over NumPy arrays and imports nothing but NumPy;
   `python_som._convert` adapts pandas at the boundary; `python_som._som` keeps the state, the
-  validation and the training loops. **No public name, signature or numerical result changes** —
+  validation and the training loops. No public name, signature or numerical result changes:
   trained weights are bit-identical to 0.3.0 for every combination of training mode and neighborhood
   function, which is asserted rather than assumed.
 
@@ -94,7 +94,7 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
   point of violation. `architecture-profile.toml` records the style.
 
 - The two update rules are now pure functions returning new arrays. An earlier note claimed this
-  was also faster; **it is not**, and the claim is withdrawn. Measured with interleaved arms and
+  was also faster. It is not, and the claim is withdrawn. Measured with interleaved arms and
   medians, the pure form is about 10% slower on small maps and indistinguishable above roughly
   100x100. The cost buys functions testable without constructing a network, and is single-digit
   milliseconds over a 10,000-iteration run on a 20x20 map. `benchmarks/bench_update.py` is the
@@ -104,19 +104,19 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 - **Batch training is 1.2x to 1.5x faster**, with identical results. Eq. (8) needs the neighborhood
   between every pair of nodes, and the previous implementation got it by evaluating the neighborhood
-  function once per node, once per iteration — 1600 full-grid evaluations per iteration on a 40x40
-  map. Profiling made that **42% of batch training**, more than the contraction it feeds.
+  function once per node, once per iteration, which is 1600 full-grid evaluations per iteration on a 40x40
+  map. Profiling made that 42% of batch training, more than the contraction it feeds.
 
   It is unnecessary, because a neighborhood depends only on the offset between two nodes and never
   on where the winner sits. One kernel over every reachable offset serves the whole grid, and
-  each node's neighborhood is a slice of it — a view, not a copy. The kernel costs `4xy` floats,
+  each node's neighborhood is a slice of it, as a view rather than a copy. The kernel costs `4xy` floats,
   111 KB on a 60x60 map, against the 800 MB a full `(x, y, x, y)` tensor would need.
 
   The offset-only dependence holds on a torus as well as a flat grid, because making the fold depend
   on the offset alone is exactly what the minimum-image convention does. It holds per axis, so mixed
   maps (one axis wrapping, one not) need no special case.
 
-  **Results do not change.** Trained weights are bit-identical for all 72 working combinations of
+  Results do not change: trained weights are bit-identical for all 72 working combinations of
   training mode, neighborhood function, initializer and cyclic setting, and the kernel is checked
   against per-node evaluation across 40,832 combinations at exactly `0.0`. The gaussian gains more
   than the bubble, whose per-node form is cheaper to begin with. `benchmarks/bench_batch.py` is the
@@ -129,16 +129,16 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 - **Linear initialization was inaccurate for data far from the origin.** Through 0.3.0 its PCA went
   to scikit-learn's `auto` solver, which since 1.5 selects `covariance_eigh` when samples outnumber
-  features — it forms the covariance matrix, and squaring the data squares its condition number. On
-  `(150, 4)` data offset by 1e7, the second explained variance was wrong by **5.8%**, and the
+  features. It forms the covariance matrix, and squaring the data squares its condition number. On
+  `(150, 4)` data offset by 1e7, the second explained variance was wrong by 5.8%, and the
   resulting models differed from the correct ones by 2.43 against a total model spread of 2.0: an
   error larger than the structure being initialized.
 
   The NumPy implementation decomposes the centred matrix directly and agrees with a `longdouble`
-  reference to ~1e-15. Data offset far from the origin is not exotic — timestamps, easting/northing
+  reference to ~1e-15. Data offset far from the origin is not exotic: timestamps, easting/northing
   coordinates and absolute sensor readings all look like it.
 
-  **This changes results.** For data near the origin the difference is floating-point noise (~1e-15,
+  This changes results. For data near the origin the difference is floating-point noise (~1e-15,
   and `auto_dimensions` is unaffected because it standardizes first). Far from the origin it is
   large, and 0.4.0 is the correct one.
 
@@ -148,11 +148,11 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 ### Removed
 
 - **pandas and scikit-learn are no longer required.** `numpy` is the only runtime dependency.
-  Measured on Linux/CPython 3.12, a fresh install drops from **333 MB across 10 packages to 69 MB
-  across 1** — a 264 MB reduction, 79% of the payload, since the two also pulled in scipy, joblib,
+  Measured on Linux/CPython 3.12, a fresh install drops from 333 MB across 10 packages to 69 MB
+  across 1, a 264 MB reduction and 79% of the payload, since the two also pulled in scipy, joblib,
   narwhals, python-dateutil, six and threadpoolctl.
 
-  Nothing is lost, and input support is *wider*. pandas was used for one `isinstance` check before
+  Nothing is lost, and input support is wider. pandas was used for one `isinstance` check before
   calling `.to_numpy()`; `np.asarray` already does that via the `__array__` protocol. Because
   that is a protocol rather than a library, polars, pyarrow, xarray and CuPy objects now work too,
   without this package knowing they exist.

README.md

@@ -5,8 +5,8 @@
 [![Python versions](https://img.shields.io/pypi/pyversions/python-som.svg)](https://pypi.org/project/python-som/)
 [![License: MIT](https://img.shields.io/badge/License-MIT-blue.svg)](https://github.com/andremsouza/python-som/blob/master/LICENSE)
 
-Implementation of Kohonen's 2-D self-organizing map, built on NumPy, with scikit-learn for PCA and
-standardization. Accepts NumPy arrays, pandas DataFrames and plain lists.
+Implementation of Kohonen's 2-D self-organizing map. NumPy is the only dependency. Accepts NumPy
+arrays, pandas DataFrames, polars, pyarrow, and anything else implementing the `__array__` protocol.
 
 [Documentation](https://andremsouza.github.io/python-som/) ·
 [Changelog](https://github.com/andremsouza/python-som/blob/master/CHANGELOG.md)
@@ -36,7 +36,7 @@ winner = som.winner(data[0])
 ```
 
 A full worked example with plots is in [examples/iris.py](https://github.com/andremsouza/python-som/blob/master/examples/iris.py) and in the
-[getting-started guide](https://andremsouza.github.io/python-som/getting-started/).
+[getting-started guide](https://andremsouza.github.io/python-som/tutorial/first-map/).
 
 ![U-matrix of a SOM trained on Iris](https://raw.githubusercontent.com/andremsouza/python-som/master/docs/assets/iris.png)
 
@@ -68,7 +68,7 @@ update of Kohonen Eq. (8) is a weighted mean whose denominator is not sign-defin
 neighborhood function. Use `mode='random'` or `mode='sequential'`; `mode='batch'` raises a
 `ValueError`.
 
-See [Neighborhood functions](https://andremsouza.github.io/python-som/neighborhood-functions/) for
+See [Neighborhood functions](https://andremsouza.github.io/python-som/reference/neighborhood-functions/) for
 the derivations, including why the Mexican hat is not an outer product of two 1-D wavelets.
 
 ## Upgrading

architecture-profile.toml

@@ -4,7 +4,7 @@
 # around it (`_som.py`, `_convert.py`) holds state, validation, progress output, and the pandas
 # adapter. See docs/ and `src/python_som/_core/__init__.py` for the rule in prose.
 #
-# The rule that matters day to day — the core imports nothing but NumPy — is enforced in CI by
+# The rule that matters day to day (the core imports nothing but NumPy) is enforced in CI by
 # ruff's TID251 (see `[tool.ruff.lint.flake8-tidy-imports.banned-api]` in pyproject.toml) and by
 # tests/test_core_boundary.py. This file is not run in CI: the checker that reads it is a
 # proprietary local review aid, not redistributable, so it cannot be a required check on a public

docs/explanation/artifact-safety.md

@@ -0,0 +1,54 @@
+# Artifact safety
+
+Loading a `.npz` written by `save_npz` cannot execute code. Loading a `.pkl` can. That is a
+difference in kind rather than degree, and it is the reason this package ships an artifact format
+at all: not to forbid `pickle`, but to make the safe path the obvious one so nobody has to reach
+for the unsafe one.
+
+What follows is the precise boundary, because "safer" without a boundary is not a claim anyone can
+act on.
+
+Safer than a pickle, and the difference is categorical rather than a matter of degree, but not
+unconditionally safe, so here is exactly what is and is not guaranteed.
+
+**Cannot happen:**
+
+- **Code execution.** `allow_pickle=False` is passed explicitly, so a crafted object array is
+  refused by NumPy, not unpickled. The test suite builds a payload that really would create a file,
+  proves it executes when NumPy is *allowed* to unpickle it, then shows the loader refusing the same
+  file with nothing created.
+- **Loading an arbitrary function by name.** Names resolve only through the registries. Nothing from
+  the file is imported, evaluated, or looked up in a module.
+- **A partly-read file being treated as valid.** The member list, the metadata JSON, the format
+  version, the weight dimensions and the weight shape against the saved config are all checked
+  before a map is constructed. Every failure raises `ArtifactError`.
+
+**Can still happen:**
+
+- **Resource exhaustion.** An `.npz` is a zip, so a hostile file can be built to decompress to
+  something enormous. Nothing here caps that.
+
+So: treat an artifact from an untrusted source the way you would a JPEG from one. It cannot run
+code, but it can still be malformed or oversized. Do not treat it as you would a signed archive.
+
+If you need integrity as well, hash the file and record the digest alongside whatever references it;
+that is outside what this format tries to do.
+
+## Why the guarantee stops there
+
+`allow_pickle=False` closes the code-execution path and nothing else. A zip archive can still be
+crafted to decompress to something enormous, and nothing in this format caps that. Capping it
+would mean inventing a size limit with no principled value, which is the kind of unsourced
+constant this package removes when it finds one.
+
+If you need integrity as well as safety, that is a different property and wants a different tool:
+hash the file and record the digest wherever the artifact is referenced.
+
+## How the claim is tested
+
+The test suite does not assert that pickles are refused. It builds a payload with a `__reduce__`
+that creates a file, proves the payload really executes when something unpickles it, and then
+shows the loader refusing a file containing it with nothing created.
+
+Proving liveness first is the part that matters. Without it the test would pass just as well
+against an inert payload, and would be evidence of nothing.

docs/training-modes.md → docs/explanation/batch-vs-stepwise.md

@@ -76,7 +76,7 @@ stepwise modes. Kohonen Section 4.1 on the radius schedule:
 
 So a sensible starting point is `neighborhood_radius` at roughly half the shorter side of the grid,
 decaying from there, with the floor described in
-[Neighborhood functions](neighborhood-functions.md#the-radius-floor) keeping it away from zero.
+[Neighborhood functions](../reference/neighborhood-functions.md#the-radius-floor) keeping it away from zero.
 
 ## Reproducibility
 

docs/explanation/why-isotropy-matters.md

@@ -0,0 +1,87 @@
+# Why isotropy matters
+
+A neighborhood function has to depend on *one* number: the distance between two nodes on the grid.
+Not on the horizontal offset and the vertical offset separately. This page is about why that
+distinction is load-bearing rather than pedantic. Getting it wrong produces a function that looks
+plausible and inverts the behaviour the mexican hat exists to provide.
+
+## What the sources require
+
+Kohonen (2013) Eq. (5) writes the neighborhood as a function of `sqdist(c, i)`, "the square of the
+geometric distance between the nodes c and i in the grid". One scalar.
+
+Vrieze (1995) is more explicit still. Figure 3, *"The 'Mexican-hat' function of lateral
+interaction"*, plots interaction against an abscissa labelled **"Lateral distance"**. The text:
+
+> When a neuron is firing, near by situated neurons are stimulated, diminishing with increasing
+> distance to the firing neuron. From a certain distance on inhibition will take place that
+> gradually vanishes when the distance becomes very large.
+
+The coefficient is written $h_{i i_c} = 1/\lVert i_c - i \rVert$, the neighborhood as
+$N_i = \{i' \mid d(i,i') \le \rho\}$, and Vrieze notes that both spaces are assumed to be metric.
+Every formulation in both papers reduces the two axis offsets to a single distance *first*, then
+applies the profile.
+
+## The mistake that looks right
+
+The gaussian is **multiplicatively separable**:
+
+$$e^{-(dx^2 + dy^2)/2\sigma^2} = e^{-dx^2/2\sigma^2} \cdot e^{-dy^2/2\sigma^2}$$
+
+So you can compute it as an outer product of two 1-D gaussians and get exactly the isotropic 2-D
+gaussian. That works, it is a normal optimization, and it tempts you into thinking the same trick
+generalizes.
+
+It does not. Separability is a property of the exponential, not of neighborhood functions. Apply the
+same construction to the Ricker wavelet and you get something else entirely:
+
+```python
+ax = (1 - dx**2 / sigma**2) * exp(-(dx**2) / (2 * sigma**2))
+ay = (1 - dy**2 / sigma**2) * exp(-(dy**2) / (2 * sigma**2))
+h = outer(ax, ay)  # not a mexican hat
+```
+
+In the diagonal quadrants both 1-D factors are negative, so their product is **positive**. The
+function stops inhibiting exactly where it is supposed to inhibit most.
+
+Measured on a 21×21 grid with $\sigma = 3$, centre at (10, 10):
+
+| quantity | outer product | correct isotropic form |
+| --- | --- | --- |
+| $h$ at $(c + 2\sigma,\, c + 2\sigma)$ | +0.165 | −0.055 |
+| $h$ at $(c + 3\sigma,\, c + 3\sigma)$ | +0.008 | −0.001 |
+| global minimum | −0.443 | −0.135 (= $-e^{-2}$, at $r = 2\sigma$) |
+| zero-crossing locus | a **cross** ($dx = \pm\sigma \cup dy = \pm\sigma$) | a **circle**, $r = \sqrt{2}\sigma$ |
+
+The separable version puts an excitatory lobe worth 16.5% of the winner's strength where the
+function should be pushing models away, and along the diagonals it never becomes negative at all.
+Its value is not a function of any metric on the grid, which is precisely the property Eq. (5)
+requires.
+
+## How the package enforces it
+
+Every neighborhood function is built from `squared_grid_distance`, which reduces the two offsets to
+one number before any profile is applied. The three shipped functions share a single implementation
+of each formula, so the per-node form and the batch kernel cannot drift apart.
+
+The test suite asserts the property directly rather than checking golden values: equal grid distance
+must give equal $h$. That assertion fails against the separable construction and passes against the
+isotropic one.
+
+## The one deliberate exception
+
+The **bubble** is a Chebyshev ball, so it is *not* isotropic under the Euclidean metric, and that is
+intentional. Vrieze's own appendix computes `b = MAX(ABS(i - w_i), ABS(j - w_j))`, a square region
+rather than a disc. Kohonen's phrasing ("up to a certain radius from the winner") reads as Euclidean,
+so the two sources genuinely differ; this package follows Vrieze and says so rather than quietly
+picking one.
+
+The consequence is worth stating because it is easy to assume otherwise: on a large enough grid,
+nodes at equal Euclidean distance can fall on opposite sides of the boundary. The smallest case is a
+radius of $\sqrt{50}$, where $(5, 5)$ lies inside a $\sigma = 5$ square and $(7, 1)$ lies outside.
+
+## Further reading
+
+- [Neighborhood functions](../reference/neighborhood-functions.md) for the formulas and constants.
+- [Batch vs stepwise](batch-vs-stepwise.md) for why a signed neighborhood cannot be used with batch
+  training at all.