From 56a8de227e612e4c6d34bb25a618f776573d1867 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Moreira=20Souza?= Date: Thu, 30 Jul 2026 20:19:43 -0300 Subject: [PATCH 1/5] test: check the neighborhood functions and training loops against MiniSom Two independent implementations of Kohonen's equations agreeing is stronger evidence than either one's own tests. That reasoning already justifies tests/test_linalg_matches_sklearn.py, which failed on first run, turned out to have the wrong reference, and then exposed a real 5.8% defect in this package's linear initialization. MiniSom is the closest peer, and the two look more interchangeable than they are: activate, winner, quantization, quantization_error and get_weights all line up while some of the mathematics behind them does not. Measured: asymptotic_decay identical, at exactly 0.0 gaussian 1.1e-16, one unit in the last place sequential training 2.8e-16 relative, after up to 200 iterations batch training 1.3e-15 relative, looser because the two accumulate Eq. (8) in a different order mexican hat differs: (1-u)e^-u here, (1-2u)e^-u there, so the zero crossing sits at sqrt(2)*sigma rather than sigma bubble differs: max(|dx|,|dy|) <= round(sigma) here, strict c-sigma < i < c+sigma there. At sigma=1, 9 nodes to 1 Both disagreements are convention choices, not defects. MiniSom's mexican hat is isotropic, so it never had the separability bug fixed here in 0.2.0. They are pinned so that neither is later "fixed" into the other. Three controls are needed before the training runs can be compared at all, and each is asserted rather than assumed. Identical injected initial models, since the two RNGs and PCA initializations differ. The sequential mode rather than the random one, because arange(T) % n and np.resize(arange(n), T) are the same sequence where the random modes are genuinely different. And a constant learning rate on train_batch_offline, which otherwise relaxes toward the Eq. (8) mean by a decaying eta rather than assigning it. Note that MiniSom's train_batch is stepwise training in sequential sample order, not Kohonen batch; train_batch_offline is the Eq. (8) implementation. minisom==2.3.6 joins the dev extra: MIT, a single minisom.py, and it declares no dependencies of its own. --- pyproject.toml | 6 + tests/test_minisom_agreement.py | 453 ++++++++++++++++++++++++++++++++ uv.lock | 8 + 3 files changed, 467 insertions(+) create mode 100644 tests/test_minisom_agreement.py diff --git a/pyproject.toml b/pyproject.toml index 0b45100..9937bed 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -41,6 +41,12 @@ dependencies = [ cli = ["tqdm>=4.66"] dev = [ "hypothesis==6.163.0", + # Not a runtime dependency. tests/test_minisom_agreement.py checks our neighborhood functions, + # decay and training loops against MiniSom's, the closest comparable implementation, the same way + # tests/test_linalg_matches_sklearn.py checks our PCA against scikit-learn's. Two independent + # implementations of Kohonen's equations agreeing is worth more than either one's own tests. + # Cheap to depend on: MIT, a single minisom.py, and it declares no dependencies of its own. + "minisom==2.3.6", "mypy==2.3.0", # Not runtime dependencies. The suite exercises the DataFrame path through the port, and # tests/test_linalg_matches_sklearn.py re-derives every PCA both ways and compares. Pinned like diff --git a/tests/test_minisom_agreement.py b/tests/test_minisom_agreement.py new file mode 100644 index 0000000..5143d9b --- /dev/null +++ b/tests/test_minisom_agreement.py @@ -0,0 +1,453 @@ +"""Check this package against MiniSom, the closest comparable implementation. + +Two independent implementations of Kohonen's equations agreeing is stronger evidence than either +one's own tests, which is the same reasoning behind ``tests/test_linalg_matches_sklearn.py``. That +file earned its place: it failed on first run, the *reference* turned out to be the wrong one, and +pursuing it exposed a real 5.8% defect in this package's linear initialization. + +The reason this file exists at all is that the two libraries look far more interchangeable than they +are. Method names line up (``activate``, ``winner``, ``quantization``, ``quantization_error``, +``get_weights``) while the mathematics behind some of them does not, so a comparison written from +the names alone would compare different algorithms and attribute the difference to something else. +What agrees and what does not, measured rather than assumed: + +=========================== ============================================================= +agrees exactly ``asymptotic_decay``, at 0.0 +agrees to round-off the gaussian (1.1e-16), sequential training (2.8e-16), + batch training (3.1e-15) +**does not agree** the mexican hat: ``(1-u)e^-u`` here, ``(1-2u)e^-u`` there, so the + zero crossing sits at sqrt(2)*sigma rather than sigma +**does not agree** the bubble: ``max(|dx|,|dy|) <= round(sigma)`` here, strict + ``c-sigma < i < c+sigma`` there. At sigma=1, 9 nodes against 1 +=========================== ============================================================= + +Both disagreements are convention choices rather than defects on either side, and MiniSom's mexican +hat is isotropic, so it does **not** have the separability bug this package fixed in 0.2.0. They are +pinned here so that nobody later "fixes" one into the other, and so the comparison published in +``docs/explanation/`` cannot quietly stop being true. + +The end-to-end tests matter most. They are what makes it legitimate to time the two libraries +against each other: if the trained models agree to round-off, a difference in wall time is +implementation and nothing else. ``benchmarks/bench_vs_minisom.py`` relies on exactly that, and this +is where the claim is rechecked on every run. + +Nothing here is timed. Timing belongs in ``benchmarks/``, where a noisy shared runner cannot turn it +into a spurious failure. + +Skipped wholesale when MiniSom is absent, so the suite still runs in an environment that installed +only the runtime dependency. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING, Any + +import numpy as np +import pytest + +import python_som +from python_som._core._decay import asymptotic_decay +from python_som._core._neighborhood import bubble, gaussian, mexican_hat + +if TYPE_CHECKING: # pragma: no cover + import numpy.typing as npt + +minisom = pytest.importorskip("minisom") + +#: Round-off scale for one evaluation of a neighborhood function. Measured at 1.1e-16, which is a +#: single unit in the last place; this leaves room for a different BLAS without admitting a +#: different formula. +NEIGHBORHOOD_TOLERANCE = 1e-14 + +#: Round-off scale for a whole training run, relative to the largest model component. Measured at +#: 1.5e-16 for sequential and 1.3e-15 for batch over 100 to 200 iterations, the latter larger +#: because the two accumulate Eq. (8) in a different order. Three orders of headroom: loose enough +#: to survive another machine, tight enough that a real divergence in the update rule fails it. +TRAINING_TOLERANCE = 1e-12 + +#: Independent of the fixture seeds: this file compares two implementations, not two runs. +SEED = 20260730 + +#: Grids to sweep, including a degenerate single-row map. +SHAPES = [(5, 5), (7, 4), (12, 9), (1, 6)] + +#: Radii to sweep. Below 1 and above the grid are both legal and both worth covering. +RADII = [0.5, 1.0, 2.5, 3.0, 7.0] + + +def _data(n_samples: int, n_features: int) -> npt.NDArray[np.floating]: + """Build a reproducible dataset. + + :param n_samples: Number of samples. + :param n_features: Number of features. + :return: The dataset. + """ + rng = np.random.default_rng(SEED) + return rng.normal(size=(n_samples, n_features)) + + +def _peer(shape: tuple[int, int], n_features: int, **kwargs: object) -> Any: # noqa: ANN401 + """Build a MiniSom of the given shape. + + ``Any`` because MiniSom ships no type information, so there is nothing more precise to say + about what comes back. + + :param shape: Grid shape. + :param n_features: Number of input features. + :param kwargs: Passed through to the MiniSom constructor. + :return: The MiniSom. + """ + return minisom.MiniSom(shape[0], shape[1], n_features, random_seed=SEED, **kwargs) + + +def _relative_difference(ours: npt.NDArray[np.floating], theirs: npt.NDArray[np.floating]) -> float: + """Return the largest absolute difference, scaled by the largest value being compared. + + Preferred over ``assert_allclose``'s ``rtol`` because a model component that happens to sit near + zero would fail a relative comparison for no reason worth failing over. + + :param ours: This package's result. + :param theirs: MiniSom's result. + :return: The scaled difference. + """ + scale = float(np.abs(theirs).max()) + difference = float(np.abs(ours - theirs).max()) + return difference / scale if scale else difference + + +# --------------------------------------------------------------------------------------------- +# The neighborhood functions +# --------------------------------------------------------------------------------------------- + + +@pytest.mark.parametrize("shape", SHAPES) +@pytest.mark.parametrize("sigma", RADII) +def test_the_gaussians_are_the_same_function(shape: tuple[int, int], sigma: float) -> None: + """Every node of every grid, not one spot check. + + The two are written differently and it is not obvious they agree. MiniSom takes the outer + product of two one-dimensional gaussians; this package exponentiates a single squared grid + distance. They coincide because the exponential is separable, + ``exp(-(dx^2+dy^2)/2s^2) == exp(-dx^2/2s^2) * exp(-dy^2/2s^2)``, which is a property of the + gaussian alone and does not carry to the other two neighborhoods. That is the whole reason the + head-to-head benchmark uses the gaussian. + """ + peer = _peer(shape, 3) + for cx in range(shape[0]): + for cy in range(shape[1]): + ours = gaussian(shape, (cx, cy), sigma, (False, False)) + theirs = peer._gaussian((cx, cy), sigma) + assert ours.shape == theirs.shape + np.testing.assert_allclose(ours, theirs, atol=NEIGHBORHOOD_TOLERANCE) + + +def test_the_mexican_hats_differ_by_a_factor_of_two_in_the_linear_term() -> None: + """A deliberate difference, pinned so it is not mistaken for a bug in either package. + + This package uses ``(1 - u) * exp(-u)`` with ``u = r^2 / 2s^2``, so the profile crosses zero at + ``r = sqrt(2) * s``. MiniSom uses ``(1 - 2u) * exp(-u)``, crossing at ``r = s``. Both are + standard ways to write a Ricker wavelet and neither is wrong; they are simply not the same + function, so the two cannot be compared at equal sigma. + + Worth stating plainly because this package's own mexican hat *was* wrong once, in the form + contributed in PR #2: an outer product of two one-dimensional Ricker wavelets, which is not + radial and produced a positive side lobe of +0.165 on the diagonal at 2 sigma where the correct + value is -0.055. MiniSom's applies its profile to a single squared distance and so never had + that defect. This test asserts the remaining difference is only the factor of two. + """ + shape, centre, sigma = (41, 41), (20, 20), 3.0 + ours = mexican_hat(shape, centre, sigma, (False, False)) + theirs = _peer(shape, 3)._mexican_hat(centre, sigma) + + radius = np.sqrt( + np.add.outer( + (np.arange(shape[0]) - centre[0]) ** 2.0, (np.arange(shape[1]) - centre[1]) ** 2.0 + ) + ) + # The crossing is bracketed rather than solved for, since the grid only samples the profile. + assert radius[ours > 0].max() < np.sqrt(2) * sigma < radius[ours < 0].min() + assert radius[theirs > 0].max() < sigma < radius[theirs < 0].min() + + # Both are radial, which is the property PR #2 got wrong, and both peak at 1 on the winner. + assert ours[centre] == pytest.approx(1.0) + assert theirs[centre] == pytest.approx(1.0) + for profile in (ours, theirs): + for distance in np.unique(radius): + shell = profile[radius == distance] + assert shell.max() - shell.min() == pytest.approx(0.0, abs=1e-15) + + +def test_the_bubbles_select_different_neighbourhoods() -> None: + """Also deliberate, and a much larger difference than the mexican hat's. + + This package rounds sigma and includes the boundary, ``max(|dx|,|dy|) <= round(sigma)``, which + at sigma=1 is the 3x3 block around the winner. MiniSom uses strict inequalities against an + unrounded sigma, ``c - sigma < i < c + sigma``, selecting only the winner at sigma=1. Nine + nodes against one is not a rounding difference, and a benchmark that used the bubble at equal + sigma would be comparing a neighbourhood against a point update. + """ + shape, centre = (9, 9), (4, 4) + # sigma=1 is an integer >= 1, so MiniSom's "sigma should be an integer" warning does not fire + # and the suite's filterwarnings=error setting is not tripped. The next test covers the case + # where it does. + peer = _peer(shape, 3, neighborhood_function="bubble", sigma=1) + + assert bubble(shape, centre, 1.0, (False, False)).sum() == 9 + assert peer._bubble(centre, 1.0).sum() == 1 + + +def test_minisom_warns_about_a_non_integer_bubble_radius() -> None: + """Recorded because this suite runs ``filterwarnings = ["error"]``. + + A future test that constructs a MiniSom with the bubble and a decayed, non-integer sigma would + fail with a ``UserWarning`` rather than an assertion, which is a confusing way to find out. This + package places no such restriction on its own bubble. + """ + with pytest.warns(UserWarning, match="sigma should be an integer"): + _peer((9, 9), 3, neighborhood_function="bubble", sigma=1.5) + + +# --------------------------------------------------------------------------------------------- +# The pieces the training loops are built from +# --------------------------------------------------------------------------------------------- + + +def test_the_asymptotic_decays_are_identical() -> None: + """Exactly ``0.0``, not a tolerance: both evaluate ``x / (1 + t / (max_t / 2))``. + + This is what lets the radius and learning-rate schedules be matched between the two libraries by + naming the same decay on each side, rather than by reimplementing one of them. + """ + peer = _peer((5, 5), 3) + for x in (0.5, 1.0, 3.0, 10.0): + for max_t in (10, 100, 1000): + for t in range(0, max_t, max(1, max_t // 7)): + assert asymptotic_decay(x, t, max_t) - peer._asymptotic_decay(x, t, max_t) == 0.0 + + +def test_the_winners_agree_for_the_same_models() -> None: + """Same best-matching unit for every sample, including the coordinate convention.""" + shape, n_features = (7, 5), 4 + data = _data(40, n_features) + weights = np.random.default_rng(SEED).normal(size=(*shape, n_features)) + + som = python_som.SOM(x=shape[0], y=shape[1], input_len=n_features, random_seed=SEED) + som._weights = weights.copy() + peer = _peer(shape, n_features) + peer._weights = weights.copy() + + for sample in data: + assert som.winner(sample) == tuple(peer.winner(sample)) + + +def test_the_winners_agree_on_a_deliberate_tie() -> None: + """Ties go to the first index in C order in both, which is arbitrary but must be the same one. + + Constructed rather than hoped for: two nodes are placed at exactly equal distance from the + sample. If the two packages broke ties differently, every later comparison would drift apart on + data with duplicate models, which is common right after random initialization from a small + dataset. + """ + shape, n_features = (3, 3), 2 + weights = np.full((*shape, n_features), 10.0) + weights[0, 2] = [1.0, 0.0] + weights[2, 0] = [1.0, 0.0] + sample = np.array([1.0, 0.0]) + + som = python_som.SOM(x=shape[0], y=shape[1], input_len=n_features, random_seed=SEED) + som._weights = weights.copy() + peer = _peer(shape, n_features) + peer._weights = weights.copy() + + assert som.winner(sample) == (0, 2) + assert tuple(peer.winner(sample)) == (0, 2) + + +def test_the_quantization_errors_agree() -> None: + """The metric the benchmark reports, so the two must define it the same way. They do. + + Both are the mean Euclidean distance from each sample to its best-matching model. Note that the + similarly named ``quantization`` does **not** agree and is not interchangeable: this package + returns one distance per sample, MiniSom returns the winning model vector per sample. Only + ``quantization_error`` is the comparable pair, which is asserted below. + """ + shape, n_features = (6, 6), 3 + data = _data(50, n_features) + weights = np.random.default_rng(SEED).normal(size=(*shape, n_features)) + + som = python_som.SOM(x=shape[0], y=shape[1], input_len=n_features, random_seed=SEED) + som._weights = weights.copy() + peer = _peer(shape, n_features) + peer._weights = weights.copy() + + assert som.quantization_error(data) == pytest.approx(peer.quantization_error(data), abs=1e-14) + assert som.quantization(data).shape == (len(data),) + assert peer.quantization(data).shape == (len(data), n_features) + + +# --------------------------------------------------------------------------------------------- +# End to end: the claim the benchmark rests on +# --------------------------------------------------------------------------------------------- + +#: Grid, sample count, feature count and iteration count per end-to-end case. +TRAINING_CASES = [((8, 8), 60, 3, 100), ((12, 9), 120, 5, 200), ((20, 20), 200, 4, 150)] + +#: Initial radius and learning rate, shared by both libraries. Chosen with the floor below in mind. +SIGMA_0 = 3.0 +LEARNING_RATE_0 = 0.5 + + +def _matched_som(shape: tuple[int, int], n_features: int) -> python_som.SOM: + """Build a map whose schedules match what MiniSom will be given. + + Everything that has to agree between the two libraries is fixed here rather than at the call + sites, so a case cannot silently opt out of one of the controls. + + :param shape: Grid shape. + :param n_features: Number of input features. + :return: The map. + """ + return python_som.SOM( + x=shape[0], + y=shape[1], + input_len=n_features, + learning_rate=LEARNING_RATE_0, + learning_rate_decay=asymptotic_decay, + neighborhood_radius=SIGMA_0, + neighborhood_radius_decay=asymptotic_decay, + neighborhood_function="gaussian", + random_seed=SEED, + ) + + +@pytest.mark.parametrize(("shape", "n_samples", "n_features", "n_iteration"), TRAINING_CASES) +def test_the_radius_floor_never_binds_in_these_cases( + shape: tuple[int, int], n_samples: int, n_features: int, n_iteration: int +) -> None: + """A precondition of the two tests below, checked rather than assumed. + + This package floors the decayed radius at ``min_neighborhood_radius`` (Kohonen Section 4.2: + the radius "should always remain, say, above half of the grid spacing"). MiniSom has no such + floor, so if it ever engaged the two would be running different schedules and the agreement + below would be measuring the wrong thing. At sigma=3 over 100 iterations the final radius is + 1.007, comfortably above the 0.5 default, but that is a property of the chosen cases rather than + of the code, so it is asserted per case. + """ + del n_samples, n_features # part of the shared case tuple, not needed here + som = _matched_som(shape, 3) + for t in range(n_iteration): + floored = som._sigma(t, n_iteration) + assert floored == asymptotic_decay(SIGMA_0, t, n_iteration), ( + f"the radius floor engaged at t={t}; pick a larger sigma or fewer iterations" + ) + + +@pytest.mark.parametrize(("shape", "n_samples", "n_features", "n_iteration"), TRAINING_CASES) +def test_sequential_training_agrees_end_to_end( + shape: tuple[int, int], n_samples: int, n_features: int, n_iteration: int +) -> None: + """Kohonen Eq. (3), driven identically on both sides, must land on the same models. + + Note which MiniSom method this is. ``train_batch`` is *stepwise* training in sequential sample + order, the counterpart of ``mode="sequential"`` here; MiniSom's Eq. (8) implementation is called + ``train_batch_offline`` and is exercised by the next test. Reaching for the name that matches + this package's would compare per-sample gradient steps against a weighted mean. + + Sequential rather than random because the index sequences then coincide exactly: + ``np.resize(np.arange(n), T)`` here against ``arange(T) % n`` there. The random modes are + genuinely different, i.i.d. draws here against a shuffle of a fixed multiset there, so no seed + reconciles them. + + Measured agreement: 2.2e-16 to 2.8e-16 relative, after up to 200 iterations. + """ + data = _data(n_samples, n_features) + initial = np.random.default_rng(SEED + 1).normal(size=(*shape, n_features)) + + som = _matched_som(shape, n_features) + som._weights = initial.copy() + som.train(data, n_iteration=n_iteration, mode="sequential") + + peer = _peer( + shape, + n_features, + sigma=SIGMA_0, + learning_rate=LEARNING_RATE_0, + decay_function="asymptotic_decay", + sigma_decay_function="asymptotic_decay", + neighborhood_function="gaussian", + ) + peer._weights = initial.copy() + peer.train_batch(data, n_iteration) + + assert _relative_difference(som.get_weights(), peer.get_weights()) < TRAINING_TOLERANCE + + +@pytest.mark.parametrize(("shape", "n_samples", "n_features", "n_iteration"), TRAINING_CASES) +def test_batch_training_agrees_end_to_end( + shape: tuple[int, int], n_samples: int, n_features: int, n_iteration: int +) -> None: + """Kohonen Eq. (8) on both sides, once MiniSom's extra learning rate is removed. + + MiniSom's ``train_batch_offline`` cites the same paper this package works from, and computes the + same neighborhood-weighted mean, but then relaxes toward it rather than assigning it: + ``w <- (1 - eta) w + eta * mean``. Eq. (8) has no step size, so the two coincide only at + ``eta == 1``. Passing ``learning_rate=1.0`` alone is not enough, because the rate is still fed + through a decay; a constant callable is what pins it. + + That difference is the reason this test exists rather than being folded into the one above. It + is exactly the kind of thing that would otherwise show up as a modest, believable quality gap in + a published comparison. + + Measured agreement: 6.0e-16 to 1.3e-15 relative, looser than the sequential case because the two + accumulate the same sum in a different order. + """ + data = _data(n_samples, n_features) + initial = np.random.default_rng(SEED + 1).normal(size=(*shape, n_features)) + + som = _matched_som(shape, n_features) + som._weights = initial.copy() + som.train(data, n_iteration=n_iteration, mode="batch") + + peer = _peer( + shape, + n_features, + sigma=SIGMA_0, + learning_rate=1.0, + decay_function=lambda rate, t, max_t: 1.0, # noqa: ARG005 + sigma_decay_function="asymptotic_decay", + neighborhood_function="gaussian", + ) + peer._weights = initial.copy() + peer.train_batch_offline(data, n_iteration) + + assert _relative_difference(som.get_weights(), peer.get_weights()) < TRAINING_TOLERANCE + + +def test_minisom_batch_without_the_pinned_rate_does_not_agree() -> None: + """The control from the previous test is load-bearing, so its absence is asserted too. + + Without it the two still look similar, which is the danger: they diverge enough to move a + published quantization error but not enough to look like a bug. Guards against someone + simplifying the constant-rate callable away because "the default is close enough". + """ + shape, n_features, n_iteration = (12, 9), 5, 100 + data = _data(120, n_features) + initial = np.random.default_rng(SEED + 1).normal(size=(*shape, n_features)) + + som = _matched_som(shape, n_features) + som._weights = initial.copy() + som.train(data, n_iteration=n_iteration, mode="batch") + + peer = _peer( + shape, + n_features, + sigma=SIGMA_0, + learning_rate=LEARNING_RATE_0, + decay_function="asymptotic_decay", + sigma_decay_function="asymptotic_decay", + neighborhood_function="gaussian", + ) + peer._weights = initial.copy() + peer.train_batch_offline(data, n_iteration) + + assert _relative_difference(som.get_weights(), peer.get_weights()) > TRAINING_TOLERANCE diff --git a/uv.lock b/uv.lock index 5d50574..2da1ee0 100644 --- a/uv.lock +++ b/uv.lock @@ -1523,6 +1523,12 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/2c/19/04f9b178c2d8a15b076c8b5140708fa6ffc5601fb6f1e975537072df5b2a/mergedeep-1.3.4-py3-none-any.whl", hash = "sha256:70775750742b25c0d8f36c55aed03d24c3384d17c951b3175d898bd778ef0307", size = 6354, upload-time = "2021-02-05T18:55:29.583Z" }, ] +[[package]] +name = "minisom" +version = "2.3.6" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/6d/9c/999bc2f1bee4cfbb6157ac03db32f5a831322479151941bc997774c405c6/minisom-2.3.6.tar.gz", hash = "sha256:1fca9cfeef62b221a89930ad3422b7a24a2bd66ae43ef97332e165f5975c09f1", size = 13140, upload-time = "2026-02-11T12:04:56.267Z" } + [[package]] name = "mkdocs" version = "1.6.1" @@ -2578,6 +2584,7 @@ cli = [ ] dev = [ { name = "hypothesis" }, + { name = "minisom" }, { name = "mypy" }, { name = "pandas", version = "2.3.3", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version < '3.11'" }, { name = "pandas", version = "3.0.5", source = { registry = "https://pypi.org/simple" }, marker = "python_full_version >= '3.11'" }, @@ -2616,6 +2623,7 @@ requires-dist = [ { name = "bandit", marker = "extra == 'analysis'", specifier = "==1.9.4" }, { name = "hypothesis", marker = "extra == 'dev'", specifier = "==6.163.0" }, { name = "matplotlib", marker = "extra == 'examples'", specifier = ">=3.8" }, + { name = "minisom", marker = "extra == 'dev'", specifier = "==2.3.6" }, { name = "mkdocs-material", marker = "extra == 'docs'", specifier = "==9.7.7" }, { name = "mkdocs-redirects", marker = "extra == 'docs'", specifier = "==1.2.2" }, { name = "mkdocstrings-python", marker = "extra == 'docs'", specifier = "==2.0.5" }, From 885ff444a9989ea7a7effbeced8b4272e44866ba Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Moreira=20Souza?= Date: Thu, 30 Jul 2026 20:36:18 -0300 Subject: [PATCH 2/5] bench: compare against MiniSom under a controlled protocol Both performance claims this package makes are against its own previous code. Neither says anything about the alternative a user is actually choosing between. The script prints two tables. "Same mathematics" applies the six controls that reduce the two libraries to the same computation, verifies the trained models agree to within 1e-12 relative, and only then reports a timing; the controls are the ones tests/test_minisom_agreement.py asserts, so the protocol cannot quietly stop holding. "Own initializer" holds the algorithm and every hyperparameter fixed and lets each library seed its own models, which isolates the one default the two do not already share. Measured on a 400-sample, 8-feature dataset, medians of 9 interleaved repeats: 20x20 batch 1.03x faster 40x40 batch 1.34x slower 60x60 batch 1.64x slower 20x20 sequential 1.12x faster 40x40 sequential 1.06x faster 60x60 sequential 1.08x faster Batch is the finding and it is not favourable. The gap tracks the ratio of nodes to samples, because batch_update walks every node in Python: 3600 nodes over 30 iterations is 108,000 einsum calls at 60x60, where MiniSom's node-side update is a single vectorised divide and its Python loop runs over the 400 samples instead. Fixing that means changing _core/_update.py, so it is not in this diff. An earlier version of this script reported sequential as 7.10x SLOWER, which was a defect in the harness rather than a result. SOM.train computes the quantization error over the whole dataset to fill in its TrainingReport and MiniSom's entry points do not, so timing train charged this package for a pass the other arm never made. On the sequential cases that pass is larger than the training itself: 30 steps touch 30 samples, the report touches all 400. Timing _train_batch and _train_stepwise directly, as bench_batch.py already does, gives the numbers above. Recorded in the module docstring because both figures were reproducible and only one of them measured training. compare() gains a repeats parameter so a cross-library case, which trains a whole map per repeat, need not use the 31 that a single update step can afford. --- benchmarks/bench_update.py | 11 +- benchmarks/bench_vs_minisom.py | 408 +++++++++++++++++++++++++++++++++ 2 files changed, 417 insertions(+), 2 deletions(-) create mode 100644 benchmarks/bench_vs_minisom.py diff --git a/benchmarks/bench_update.py b/benchmarks/bench_update.py index 104cfdc..4f37713 100644 --- a/benchmarks/bench_update.py +++ b/benchmarks/bench_update.py @@ -117,16 +117,23 @@ def run_pure(case: Case) -> npt.NDArray[np.floating]: return weights -def compare(first: Callable[[], object], second: Callable[[], object]) -> tuple[float, float, bool]: +def compare( + first: Callable[[], object], second: Callable[[], object], repeats: int = REPEATS +) -> tuple[float, float, bool]: """Time two callables with their repeats interleaved. + Also used by ``bench_vs_minisom.py``, which is why ``repeats`` is a parameter: a cross-library + case trains a whole map per repeat, so it cannot afford the 31 this script uses for a single + update step. Everything else about the method is deliberately shared rather than reimplemented. + :param first: One arm. :param second: The other arm. + :param repeats: Number of interleaved repeats. Keep it odd, so the median is an observation. :return: Median of each, and whether their interquartile ranges overlap. """ first_times: list[float] = [] second_times: list[float] = [] - for _ in range(REPEATS): + for _ in range(repeats): first_times.append(timeit.timeit(first, number=1)) second_times.append(timeit.timeit(second, number=1)) diff --git a/benchmarks/bench_vs_minisom.py b/benchmarks/bench_vs_minisom.py new file mode 100644 index 0000000..9ab2b16 --- /dev/null +++ b/benchmarks/bench_vs_minisom.py @@ -0,0 +1,408 @@ +"""Compare this package against MiniSom, the closest comparable implementation. + +Run it directly; it is not part of the test suite, because a timing assertion on shared CI hardware +would be flaky:: + + uv run python benchmarks/bench_vs_minisom.py + +Method is the same as ``bench_update.py`` and ``bench_batch.py``, and the interleaving helper is +imported from the first rather than copied: **interleaved** arms so thermal and load drift is split +evenly rather than attributed to one of them, **medians with an interquartile range** rather than +minima, and **equality asserted first**. + +That last rule is doing much more work here than in the other two scripts. Within one package it +guards against a refactor that changed an answer; across two packages it is the entire basis for the +comparison being meaningful at all. MiniSom and this package look far more interchangeable than they +are, so the numbers are only worth printing once the two have been driven into computing the same +thing. Six controls do that, and every one of them is also asserted in +``tests/test_minisom_agreement.py`` so the protocol cannot quietly stop holding: + +1. The same initial models are injected into both. Seeding cannot do it: the two use different + generators and their PCA initializations are different functions. +2. The gaussian only. It is the one neighborhood the two define identically; their mexican hats + differ by a factor of two in the linear term and their bubbles select different node sets. +3. The same radius schedule, ``asymptotic_decay``, which is the same formula in both packages. +4. Cases chosen so this package's radius floor never engages, since MiniSom has no such floor. + Asserted, because it is a property of the chosen numbers rather than of the code. +5. For batch, MiniSom's learning rate pinned to a constant 1.0. ``train_batch_offline`` relaxes + toward the Eq. (8) mean by a decaying ``eta`` where Eq. (8) simply assigns it. +6. Sequential rather than random sample order, because ``arange(T) % n`` and + ``np.resize(arange(n), T)`` are the same sequence. The random modes are genuinely different + (i.i.d. draws here, a shuffled fixed multiset there) and no seed reconciles them. + +Note which MiniSom methods those map to. ``train_batch`` is *stepwise* training in sequential sample +order, the counterpart of ``mode="sequential"`` here. Kohonen Eq. (8) is ``train_batch_offline``. +Reaching for the name that matches this package's would time a weighted mean against per-sample +gradient steps. + +**What is timed.** The training loops, ``_train_batch`` and ``_train_stepwise``, not the public +``train``. This is not a shortcut, it is the difference between a fair comparison and a wrong one. +``train`` also computes the quantization error over the whole dataset to fill in its +``TrainingReport``; MiniSom's entry points do nothing of the kind. Timing ``train`` therefore +charges this package for a full extra pass over the data that the other arm never pays for, and on +the +sequential cases that pass is larger than the training itself: 30 steps touch 30 samples, the report +touches all 400. The first version of this script did exactly that and reported this package as +**7.10x slower** on the largest sequential case. Timing the loop, the same case is **1.08x faster**. +Both numbers were reproducible; only one of them measured training. ``bench_batch.py`` avoids the +same trap by reproducing the loop rather than calling the public method. + +Two tables are printed and the difference between them matters: + +- **Same mathematics** applies all six controls, so both libraries compute the same result and the + only thing left to measure is how fast they compute it. +- **Own initializer** keeps the algorithm and every hyperparameter fixed and lets each library seed + its own models the way its documentation says to. It is not an algorithm comparison; it measures + one packaging decision, which happens to be the only default the two do not already share. +""" + +from __future__ import annotations + +import functools +import tracemalloc +from importlib.metadata import version +from typing import TYPE_CHECKING + +import numpy as np +from bench_update import compare +from minisom import MiniSom + +import python_som +from python_som._core._decay import asymptotic_decay +from python_som._core._distance import euclidean_distance +from python_som._core._match import quantization + +if TYPE_CHECKING: # pragma: no cover + from collections.abc import Callable + + import numpy.typing as npt + +#: Repeats per arm. Odd, so the median is an observation rather than an average of two. +REPEATS = 9 + +#: Grid, sample count, feature count and iteration count per case. +CASES = [((20, 20), 200, 4, 60), ((40, 40), 300, 6, 40), ((60, 60), 400, 8, 30)] + +#: Initial radius, shared by both libraries. Large enough that control 4 holds for every case above. +SIGMA_0 = 3.0 + +#: Initial learning rate, shared by both libraries. Unused by batch training. +LEARNING_RATE_0 = 0.5 + +#: Largest tolerated disagreement between the trained models, relative to the largest component. +#: Measured at 2.8e-16 for sequential and 1.3e-15 for batch, so this is three orders of headroom. +#: Exact equality is not available across two packages: they accumulate the same sums in different +#: orders. It is the same constant as ``tests/test_minisom_agreement.py``, for the same reason. +TOLERANCE = 1e-12 + +#: Fixed so the reported numbers can be reproduced. +SEED = 20260730 + + +def build( + shape: tuple[int, int], n_samples: int, n_features: int +) -> tuple[npt.NDArray[np.floating], npt.NDArray[np.floating]]: + """Build the dataset and the shared initial models for one case. + + :param shape: Grid shape. + :param n_samples: Number of samples. + :param n_features: Number of input features. + :return: The dataset and the initial models. + """ + rng = np.random.default_rng(SEED) + return rng.normal(size=(n_samples, n_features)), rng.normal(size=(*shape, n_features)) + + +def our_error(weights: npt.NDArray[np.floating], data: npt.NDArray[np.floating]) -> float: + """Return the quantization error of any library's models, under this package's definition. + + Applied to every arm rather than calling each library's own method. Here the two definitions do + agree, which ``tests/test_minisom_agreement.py`` asserts, but recomputing keeps the reported + column meaning one thing regardless of whose models produced it. + + :param weights: Models, of shape ``(x, y, n_features)``. + :param data: Dataset of shape ``(n_samples, n_features)``. + :return: Mean distance from each sample to its best-matching model. + """ + return float(quantization(data, weights, euclidean_distance).mean()) + + +# --------------------------------------------------------------------------------------------- +# The two arms +# --------------------------------------------------------------------------------------------- + + +def train_ours( + shape: tuple[int, int], + data: npt.NDArray[np.floating], + initial: npt.NDArray[np.floating], + n_iteration: int, + mode: str, +) -> npt.NDArray[np.floating]: + """Train this package's map from the injected models. + + :param shape: Grid shape. + :param data: Training dataset. + :param initial: Models to start from. + :param n_iteration: Number of iterations. + :param mode: Either ``'batch'`` or ``'sequential'``. + :return: The trained models. + """ + som = python_som.SOM( + x=shape[0], + y=shape[1], + input_len=data.shape[1], + learning_rate=LEARNING_RATE_0, + learning_rate_decay=asymptotic_decay, + neighborhood_radius=SIGMA_0, + neighborhood_radius_decay=asymptotic_decay, + neighborhood_function="gaussian", + random_seed=SEED, + ) + som._weights = initial.copy() # noqa: SLF001 control 1: the same models on both sides + # The private loops rather than `train`. See "What is timed" in the module docstring. + if mode == "batch": + som._train_batch(data, n_iteration, verbose=False) # noqa: SLF001 + else: + som._train_stepwise(data, n_iteration, mode=mode, verbose=False) # noqa: SLF001 + return som.get_weights() + + +def train_theirs( + shape: tuple[int, int], + data: npt.NDArray[np.floating], + initial: npt.NDArray[np.floating], + n_iteration: int, + mode: str, +) -> npt.NDArray[np.floating]: + """Train MiniSom from the injected models, with the schedules matched. + + :param shape: Grid shape. + :param data: Training dataset. + :param initial: Models to start from. + :param n_iteration: Number of iterations. + :param mode: Either ``'batch'`` or ``'sequential'``. + :return: The trained models. + """ + batch = mode == "batch" + peer = MiniSom( + shape[0], + shape[1], + data.shape[1], + sigma=SIGMA_0, + # Control 5: Eq. (8) has no step size, so the blend has to be turned off for batch. For + # stepwise the rate is part of Eq. (3) and matches this package's default. + learning_rate=1.0 if batch else LEARNING_RATE_0, + decay_function=(lambda rate, t, max_t: 1.0) if batch else "asymptotic_decay", # noqa: ARG005 + sigma_decay_function="asymptotic_decay", + neighborhood_function="gaussian", + random_seed=SEED, + ) + peer._weights = initial.copy() # noqa: SLF001 + if batch: + peer.train_batch_offline(data, n_iteration) + else: + peer.train_batch(data, n_iteration) # despite the name, this is stepwise sequential + weights: npt.NDArray[np.floating] = peer.get_weights() + return weights + + +def seed_and_train_ours( + shape: tuple[int, int], + data: npt.NDArray[np.floating], + n_iteration: int, +) -> npt.NDArray[np.floating]: + """Seed and train with this package's own initializer and defaults. + + :param shape: Grid shape. + :param data: Training dataset. + :param n_iteration: Number of iterations. + :return: The trained models. + """ + som = python_som.SOM(x=shape[0], y=shape[1], input_len=data.shape[1], random_seed=SEED) + som.weight_initialization(mode="linear", data=data) + som._train_batch(data, n_iteration, verbose=False) # noqa: SLF001 see "What is timed" + return som.get_weights() + + +def seed_and_train_theirs( + shape: tuple[int, int], + data: npt.NDArray[np.floating], + n_iteration: int, +) -> npt.NDArray[np.floating]: + """Seed and train with MiniSom's own initializer and defaults. + + :param shape: Grid shape. + :param data: Training dataset. + :param n_iteration: Number of iterations. + :return: The trained models. + """ + peer = MiniSom(shape[0], shape[1], data.shape[1], random_seed=SEED) + peer.pca_weights_init(data) + peer.train_batch_offline(data, n_iteration) + weights: npt.NDArray[np.floating] = peer.get_weights() + return weights + + +def assert_the_floor_never_engaged(shape: tuple[int, int], n_iteration: int) -> None: + """Check control 4 for one case. + + This package floors the decayed radius at ``min_neighborhood_radius``; MiniSom does not. If the + floor engaged, the two would be running different schedules and the agreement check below would + pass or fail for the wrong reason. + + :param shape: Grid shape. + :param n_iteration: Number of iterations. + :raises AssertionError: If the floor changed the radius at any step. + """ + som = python_som.SOM( + x=shape[0], y=shape[1], input_len=2, neighborhood_radius=SIGMA_0, random_seed=SEED + ) + for t in range(n_iteration): + if som._sigma(t, n_iteration) != asymptotic_decay(SIGMA_0, t, n_iteration): # noqa: SLF001 + message = f"{shape}: the radius floor engaged at t={t}, so the schedules differ" + raise AssertionError(message) + + +def peak_megabytes(run: Callable[[], object]) -> float: + """Return the peak memory one call allocates, in megabytes. + + Measured in a pass of its own rather than during the timed repeats, because tracing every + allocation slows the thing being timed. NumPy registers its allocator with ``tracemalloc``, so + array data is included and not just the Python objects wrapping it. + + :param run: The callable to measure. + :return: Peak traced allocation in megabytes. + """ + tracemalloc.start() + try: + run() + return tracemalloc.get_traced_memory()[1] / 1e6 + finally: + tracemalloc.stop() + + +# --------------------------------------------------------------------------------------------- +# The tables +# --------------------------------------------------------------------------------------------- + + +def same_mathematics() -> list[str]: + """Time both libraries on identical inputs, after proving they compute the same thing. + + :return: The rendered table. + :raises AssertionError: If the two disagree by more than the tolerance. + """ + header = ( + f"{'map':>9} {'samples':>8} {'features':>9} {'iters':>6} {'mode':>11} " + f"{'python-som':>11} {'MiniSom':>11} {'ratio':>8} {'peak MB':>16} verdict" + ) + lines = [header, "-" * len(header)] + + for shape, n_samples, n_features, n_iteration in CASES: + data, initial = build(shape, n_samples, n_features) + assert_the_floor_never_engaged(shape, n_iteration) + + for mode in ("batch", "sequential"): + ours = functools.partial(train_ours, shape, data, initial, n_iteration, mode) + theirs = functools.partial(train_theirs, shape, data, initial, n_iteration, mode) + + trained_ours, trained_theirs = ours(), theirs() + difference = float(np.abs(trained_ours - trained_theirs).max()) / float( + np.abs(trained_theirs).max() + ) + if difference > TOLERANCE: + message = ( + f"{shape} {mode}: the two libraries disagree by {difference:.2e} relative, " + f"above the {TOLERANCE:.0e} tolerance, so timing them proves nothing" + ) + raise AssertionError(message) + + peak_ours, peak_theirs = peak_megabytes(ours), peak_megabytes(theirs) + median_ours, median_theirs, overlap = compare(ours, theirs, REPEATS) + ratio = median_theirs / median_ours + if overlap: + verdict = "indistinguishable (IQRs overlap)" + elif ratio > 1: + verdict = f"python-som {ratio:.2f}x faster" + else: + verdict = f"python-som {1 / ratio:.2f}x slower" + + lines.append( + f"{shape[0]:>4}x{shape[1]:<4} {n_samples:>8} {n_features:>9} {n_iteration:>6} " + f"{mode:>11} {median_ours * 1e3:>9.1f}ms {median_theirs * 1e3:>9.1f}ms " + f"{ratio:>7.2f}x {peak_ours:>7.1f} /{peak_theirs:>7.1f} {verdict}" + ) + + lines.append( + f"\nmedian of {REPEATS} interleaved repeats; trained models verified equal to within " + f"{TOLERANCE:.0e} relative before timing. Peak MB is python-som / MiniSom, measured " + f"separately from the timings." + ) + return lines + + +def own_initializer() -> list[str]: + """Compare the two initializers, with the algorithm and every hyperparameter held fixed. + + The defaults happen to be identical on both sides (sigma 1.0, learning rate 0.5, asymptotic + decay on both schedules, gaussian, euclidean), so the initializer is the only default that + differs and this table isolates it rather than confounding it with anything else. + + The two are genuinely different functions, not two spellings of one. This package places models + at ``mean + c1*sqrt(lambda1)*v1 + c2*sqrt(lambda2)*v2``, scaling each principal direction by the + data's actual extent along it, from an SVD of the centred matrix. MiniSom places them at + ``mean + c1*v1 + c2*v2`` with unit eigenvectors, so the initial sheet is one unit across + whatever the data's spread, and it eigendecomposes the covariance matrix rather than the data. + + :return: The rendered table. + """ + header = ( + f"{'map':>9} {'samples':>8} {'features':>9} {'iters':>6} " + f"{'python-som':>11} {'MiniSom':>11} {'qe ours':>9} {'qe theirs':>10} better" + ) + lines = [header, "-" * len(header)] + + for shape, n_samples, n_features, n_iteration in CASES: + rng = np.random.default_rng(SEED) + data = rng.normal(size=(n_samples, n_features)) * 10.0 + 100.0 + + # functools.partial rather than a closure: a closure over the loop variables would time + # whatever they held when it ran, not when it was written. Same reason bench_batch.py does. + ours = functools.partial(seed_and_train_ours, shape, data, n_iteration) + theirs = functools.partial(seed_and_train_theirs, shape, data, n_iteration) + + error_ours = our_error(ours(), data) + error_theirs = our_error(theirs(), data) + median_ours, median_theirs, _ = compare(ours, theirs, REPEATS) + + better = "python-som" if error_ours < error_theirs else "MiniSom" + lines.append( + f"{shape[0]:>4}x{shape[1]:<4} {n_samples:>8} {n_features:>9} {n_iteration:>6} " + f"{median_ours * 1e3:>9.1f}ms {median_theirs * 1e3:>9.1f}ms " + f"{error_ours:>9.4f} {error_theirs:>10.4f} {better}" + ) + + lines.append( + "\nBatch training on both sides; only the initializer differs. Quantization error is " + "recomputed under this package's definition for both arms. Data is deliberately offset " + "from the origin and scaled by 10, which is where the two initializers diverge most." + ) + return lines + + +def main() -> None: + """Print both tables, with the versions they were measured against.""" + print( + f"python-som {python_som.__version__}, MiniSom {version('minisom')}, " + f"NumPy {np.__version__}\n" + ) + print("SAME MATHEMATICS: both libraries computing the same result, timed") + print("\n".join(same_mathematics())) + print("\n") + print("OWN INITIALIZER: same algorithm and hyperparameters, each library seeding its own") + print("\n".join(own_initializer())) + + +if __name__ == "__main__": + main() From 020c0be6d483126b2775a3fa6bc8619d9613b1c0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Moreira=20Souza?= Date: Thu, 30 Jul 2026 22:36:48 -0300 Subject: [PATCH 3/5] bench: compare against SOMPY in an isolated environment SOMPY's batch training is Kohonen Eq. (8) by way of the Helsinki SOM Toolbox, the same algorithm this package implements, so comparing the two is comparing implementations of one published equation. Getting there costs more machinery than MiniSom did, for one reason: SOMPY cannot be installed alongside this package at all. It uses np.Inf at four sites, three of them as default arguments evaluated when the class body executes. np.Inf was removed in NumPy 2.0, so `import sompy` raises AttributeError on every NumPy this package supports. Two further reasons to keep it at arm's length regardless: importing it reconfigures the root logger for the whole process, and it pulls in matplotlib and scikit-image through sompy.visualization. So the comparison runs across a process boundary. sompy_worker.py trains one map inside a hand-built .venv-sompy and times only that call; bench_vs_sompy.py launches one worker per repeat, which keeps the arms interleaved while leaving the half second of interpreter start-up outside every measurement. If the environment is absent the script prints the setup command and exits 0. Nothing installs anything. Four things about that environment were each found by an install failing, and all four are recorded next to the command: NumPy must be pinned below 2 in the same resolve as everything else, since SOMPY's own metadata says numpy >= 1.7 and anything installed afterwards drags NumPy 2 back in; Python 3.12 is the ceiling for NumPy 1.26 wheels; scikit-image, joblib and matplotlib are imported but absent from install_requires; and the pin is a commit SHA because the repository has no tags and is not on PyPI under a name that resolves to it. Measured, batch and gaussian only, medians of 9 interleaved repeats: 20x20 1.30x faster 40x40 1.74x faster 60x60 1.99x faster Trained models agree to 1.4e-07 relative and the quantization errors match to four decimals. That is looser than the 1e-12 MiniSom is held to, because SOMPY rounds its codebook to six decimals after every update and expands ||x-w||^2 as -2x.w + ||w||^2, the less stable arrangement. The measured figure is printed in the table rather than only asserted, so a reader sees what was tolerated. Two things the harness must not do, both recorded in the code. SOMPY's calculate_quantization_error returns mean(|x - m|) over every feature, an elementwise MAE rather than the mean Euclidean distance per sample, so the metric is recomputed here for both arms; SOMPY is inconsistent with itself about this, since the value it logs per epoch is an L2 distance. And the worker reports its own NumPy version in its output, because this process cannot read it: importlib.metadata would name the version SOMPY cannot run on. Cases stop at 60x60. SOMPY materialises the full (nnodes, nnodes) neighborhood, 104 MB there and 800 MB at 100x100, which is the tensor this package rejected in 0.3.0. The cap and its reason are printed with the results rather than left silent. --- .gitignore | 3 + benchmarks/bench_vs_sompy.py | 308 +++++++++++++++++++++++++++++++++++ benchmarks/sompy_worker.py | 124 ++++++++++++++ 3 files changed, 435 insertions(+) create mode 100644 benchmarks/bench_vs_sompy.py create mode 100644 benchmarks/sompy_worker.py diff --git a/.gitignore b/.gitignore index 5674cdd..178b3f2 100644 --- a/.gitignore +++ b/.gitignore @@ -109,6 +109,9 @@ venv/ ENV/ env.bak/ venv.bak/ +# Built by hand for benchmarks/bench_vs_sompy.py. SOMPY cannot be installed alongside this package, +# because it uses np.Inf and so cannot be imported on NumPy 2, so it gets an interpreter of its own. +.venv-sompy/ # Spyder project settings .spyderproject diff --git a/benchmarks/bench_vs_sompy.py b/benchmarks/bench_vs_sompy.py new file mode 100644 index 0000000..179facb --- /dev/null +++ b/benchmarks/bench_vs_sompy.py @@ -0,0 +1,308 @@ +r"""Compare this package against SOMPY, across a process boundary because it has to be. + +Run it directly, once the environment below exists:: + + uv run python benchmarks/bench_vs_sompy.py + +**SOMPY needs an environment of its own, and you have to build it by hand.** Nothing here installs +anything; if the environment is missing the script prints this command and stops:: + + uv venv .venv-sompy --python 3.12 + uv pip install --python .venv-sompy \\ + "numpy==1.26.4" "scipy==1.13.1" "scikit-learn==1.5.2" "scikit-image==0.24.0" \\ + "numexpr==2.10.1" "joblib==1.4.2" "matplotlib==3.9.2" \\ + "SOMPY @ git+https://github.com/sevamoo/SOMPY@6aca604b06e5eea1391ecf507810c7aabafc3f8b" + +Every part of that is load-bearing, and each was found by the install failing: + +- **NumPy below 2.** SOMPY uses ``np.Inf``, removed in NumPy 2.0, at class-definition time. Its own + metadata says ``numpy >= 1.7`` with no upper bound, so a plain install resolves a NumPy it cannot + import. Everything else is pinned alongside it in one resolve, because installing any of them + afterwards pulls NumPy 2 straight back in. +- **Python 3.12.** The ceiling for NumPy 1.26 wheels. +- **scikit-image, joblib and matplotlib.** Imported by SOMPY, absent from its ``install_requires``. + ``sompy/__init__`` reaches all three through ``from .visualization import *``. +- **A commit SHA.** The repository has no tags, and it is not on PyPI under a name that resolves to + it: ``pip install sompy`` gets an unrelated 2016 package by a different author. + +Which is a lot of scaffolding, so it is worth saying what it buys: SOMPY's batch training is +Kohonen Eq. (8) by way of the Helsinki SOM Toolbox, the same algorithm this package implements, and +comparing two implementations of one published equation is the most informative benchmark available. + +Method matches the other scripts in this directory: interleaved arms, medians, and **agreement +asserted before any timing is printed**. Interleaving survives the process boundary because a worker +is launched per repeat and times only its own training call, so the roughly half-second of +interpreter start-up, scipy import and matplotlib import falls outside every measurement. + +Controls, beyond the ones in ``bench_vs_minisom.py``: + +- SOMPY's normalization is turned off. At its ``'var'`` default it would train on z-scored data. +- The radius follows ``linspace(start, end, n_iteration)`` on both sides, because that is the only + schedule ``_batchtrain`` can run and this package accepts any callable. +- Batch only. SOMPY implements no stepwise training. + +**Agreement here is looser than with MiniSom, for a reason in SOMPY's code rather than in ours.** It +rounds the codebook to six decimals after every update, and its best-matching-unit search expands +``||x-w||^2`` as ``-2x.w + ||w||^2``, which is the less numerically stable arrangement. So the +tolerance is set from measurement rather than from round-off, and the measured figure is printed +with the table so a reader can see what was tolerated rather than take the word for it. + +Case sizes stop at 60x60 because SOMPY materialises the full ``(nnodes, nnodes)`` neighborhood: 104 +MB there, and 800 MB at 100x100. That cap is printed with the results rather than left silent. +""" + +from __future__ import annotations + +import json +import shutil +import statistics +import subprocess +import tempfile +import timeit +from pathlib import Path +from typing import TYPE_CHECKING + +import numpy as np + +import python_som +from python_som._core._distance import euclidean_distance +from python_som._core._match import quantization + +if TYPE_CHECKING: # pragma: no cover + import numpy.typing as npt + +#: The interpreter that can import SOMPY. Never created by this script. +SOMPY_PYTHON = Path(__file__).resolve().parent.parent / ".venv-sompy" / "bin" / "python" + +#: The worker that runs inside it. +WORKER = Path(__file__).resolve().parent / "sompy_worker.py" + +#: Printed when the environment is missing, rather than installing anything. +SETUP = """\ +SOMPY needs an environment of its own, because it cannot be imported on NumPy 2. +Build it, then run this script again: + + uv venv .venv-sompy --python 3.12 + uv pip install --python .venv-sompy \\ + "numpy==1.26.4" "scipy==1.13.1" "scikit-learn==1.5.2" "scikit-image==0.24.0" \\ + "numexpr==2.10.1" "joblib==1.4.2" "matplotlib==3.9.2" \\ + "SOMPY @ git+https://github.com/sevamoo/SOMPY@6aca604b06e5eea1391ecf507810c7aabafc3f8b" +""" + +#: Repeats per arm. Odd, so the median is an observation rather than an average of two. +REPEATS = 9 + +#: Grid, sample count, feature count and iteration count per case. Capped at 60x60; see the module +#: docstring for why that is SOMPY's limit rather than an arbitrary choice. +CASES = [((20, 20), 200, 4, 40), ((40, 40), 300, 6, 30), ((60, 60), 400, 8, 20)] + +#: The radius ramp, shared by both arms. ``_batchtrain`` runs ``linspace(start, end, n_iteration)`` +#: and nothing else, so this package is the one that has to match. The end sits above the default +#: 0.5 radius floor, so that floor never engages. +RADIUS_START = 3.0 +RADIUS_END = 1.0 + +#: Largest tolerated disagreement, relative to the largest model component. Measured at 1.4e-07 to +#: 1.9e-07 across the cases above, consistent with SOMPY rounding its codebook to six decimals on +#: every update: models of order 1 lose about 5e-07 that way. Two orders of headroom. This cannot be +#: the 1e-12 that MiniSom is held to, and the reason is in SOMPY's code rather than in this package. +TOLERANCE = 1e-5 + +#: Fixed so the reported numbers can be reproduced. +SEED = 20260730 + + +def ramp(x: float, t: int, max_t: int) -> float: + """Return the radius at iteration ``t``, matching ``np.linspace(start, end, max_t)[t]``. + + :param x: Ignored. The ramp is defined by its endpoints, not by a starting value. + :param t: Current iteration. + :param max_t: Total number of iterations. + :return: The radius for this iteration. + """ + del x + if max_t == 1: + return RADIUS_START + return RADIUS_START + (RADIUS_END - RADIUS_START) * t / (max_t - 1) + + +def our_error(weights: npt.NDArray[np.floating], data: npt.NDArray[np.floating]) -> float: + """Return the quantization error of any library's models, under this package's definition. + + Recomputed rather than read from SOMPY, which is not optional here as it was with MiniSom: + ``calculate_quantization_error`` returns ``mean(|x - m|)`` over every feature, an elementwise + mean absolute error, where the standard definition is the mean Euclidean distance per sample. + SOMPY is inconsistent with itself about this, since the value it logs per epoch during training + *is* an L2 distance. + + :param weights: Models, of shape ``(x, y, n_features)``. + :param data: Dataset of shape ``(n_samples, n_features)``. + :return: Mean distance from each sample to its best-matching model. + """ + return float(quantization(data, weights, euclidean_distance).mean()) + + +def train_ours( + shape: tuple[int, int], + data: npt.NDArray[np.floating], + initial: npt.NDArray[np.floating], + n_iteration: int, +) -> npt.NDArray[np.floating]: + """Train this package's map from the injected models, on the matched radius ramp. + + ``_train_batch`` rather than ``train``, for the reason set out in ``bench_vs_minisom.py``: + ``train`` also scores the whole dataset to fill in its report, and the other arm does not. + + :param shape: Grid shape. + :param data: Training dataset. + :param initial: Models to start from. + :param n_iteration: Number of iterations. + :return: The trained models. + """ + som = python_som.SOM( + x=shape[0], + y=shape[1], + input_len=data.shape[1], + neighborhood_radius=RADIUS_START, + neighborhood_radius_decay=ramp, + neighborhood_function="gaussian", + random_seed=SEED, + ) + som._weights = initial.copy() # noqa: SLF001 + som._train_batch(data, n_iteration, verbose=False) # noqa: SLF001 + return som.get_weights() + + +def train_theirs( + shape: tuple[int, int], + data: npt.NDArray[np.floating], + initial: npt.NDArray[np.floating], + n_iteration: int, +) -> tuple[npt.NDArray[np.floating], float, str]: + """Run one SOMPY training call in its own interpreter and bring the result back. + + The worker times its own training call, so what this function costs in interpreter start-up and + imports, roughly half a second, is not part of the number it returns. + + :param shape: Grid shape. + :param data: Training dataset. + :param initial: Models to start from. + :param n_iteration: Number of iterations. + :return: The trained models, the seconds the worker spent training, and its NumPy version. + :raises RuntimeError: If the worker fails. + """ + workdir = Path(tempfile.mkdtemp(prefix="sompy-bench-")) + try: + (workdir / "spec.json").write_text( + json.dumps( + { + "shape": list(shape), + "n_iteration": n_iteration, + "radius_start": RADIUS_START, + "radius_end": RADIUS_END, + } + ), + encoding="utf-8", + ) + np.savez(workdir / "input.npz", data=data, initial=initial) + + finished = subprocess.run( # noqa: S603 + [str(SOMPY_PYTHON), str(WORKER), str(workdir)], + capture_output=True, + text=True, + check=False, + ) + if finished.returncode != 0: + message = f"the SOMPY worker failed:\n{finished.stdout}\n{finished.stderr}" + raise RuntimeError(message) + + with np.load(workdir / "output.npz") as payload: + return payload["weights"], float(payload["seconds"]), str(payload["numpy_version"]) + finally: + shutil.rmtree(workdir, ignore_errors=True) + + +def main() -> None: + """Measure every case and print the comparison, or explain how to enable it.""" + if not SOMPY_PYTHON.exists(): + print(SETUP) + return + + header = ( + f"{'map':>9} {'samples':>8} {'features':>9} {'iters':>6} " + f"{'python-som':>11} {'SOMPY':>11} {'ratio':>8} {'max diff':>10} " + f"{'qe ours':>9} {'qe theirs':>10} verdict" + ) + lines = [header, "-" * len(header)] + # Filled in by the first worker run. Printed at the end rather than the start because this + # process cannot see the other environment's NumPy: `importlib.metadata` would report the one + # imported here, which is exactly the version SOMPY cannot run on. + their_numpy = "unknown" + + for shape, n_samples, n_features, n_iteration in CASES: + rng = np.random.default_rng(SEED) + data = rng.normal(size=(n_samples, n_features)) + initial = rng.normal(size=(*shape, n_features)) + + ours = train_ours(shape, data, initial, n_iteration) + theirs, _, their_numpy = train_theirs(shape, data, initial, n_iteration) + difference = float(np.abs(ours - theirs).max()) / float(np.abs(theirs).max()) + if difference > TOLERANCE: + message = ( + f"{shape}: the two libraries disagree by {difference:.2e} relative, above the " + f"{TOLERANCE:.0e} tolerance, so timing them proves nothing" + ) + raise AssertionError(message) + + # Interleaved across the process boundary: one worker launch per repeat, alternating with + # this package's arm, so drift is split evenly rather than attributed to one of them. + our_times: list[float] = [] + their_times: list[float] = [] + for _ in range(REPEATS): + our_times.append( + timeit.timeit(lambda: train_ours(shape, data, initial, n_iteration), number=1) # noqa: B023 + ) + their_times.append(train_theirs(shape, data, initial, n_iteration)[1]) + + median_ours = statistics.median(our_times) + median_theirs = statistics.median(their_times) + low_ours, high_ours = np.percentile(our_times, [25, 75]) + low_theirs, high_theirs = np.percentile(their_times, [25, 75]) + overlap = not (high_ours < low_theirs or high_theirs < low_ours) + + ratio = median_theirs / median_ours + if overlap: + verdict = "indistinguishable (IQRs overlap)" + elif ratio > 1: + verdict = f"python-som {ratio:.2f}x faster" + else: + verdict = f"python-som {1 / ratio:.2f}x slower" + + lines.append( + f"{shape[0]:>4}x{shape[1]:<4} {n_samples:>8} {n_features:>9} {n_iteration:>6} " + f"{median_ours * 1e3:>9.1f}ms {median_theirs * 1e3:>9.1f}ms {ratio:>7.2f}x " + f"{difference:>10.1e} {our_error(ours, data):>9.4f} " + f"{our_error(theirs, data):>10.4f} {verdict}" + ) + + nodes = CASES[-1][0][0] * CASES[-1][0][1] + lines.insert( + 0, + f"python-som {python_som.__version__} on NumPy {np.__version__}, against " + f"SOMPY @6aca604 on NumPy {their_numpy}\n", + ) + lines.append( + f"\nmedian of {REPEATS} interleaved repeats, batch training only, gaussian only. " + f"'max diff' is the agreement actually measured, against a {TOLERANCE:.0e} tolerance; " + f"SOMPY rounds its codebook to 6 decimals every iteration, so exact agreement is not " + f"available. Quantization error is recomputed under this package's definition for both " + f"arms, because SOMPY's own method returns an elementwise MAE instead.\n" + f"Capped at {CASES[-1][0][0]}x{CASES[-1][0][1]}: SOMPY builds the full " + f"({nodes}, {nodes}) neighborhood matrix, {nodes * nodes * 8 / 1e6:.0f} MB here and " + f"800 MB at 100x100." + ) + print("\n".join(lines)) + + +if __name__ == "__main__": + main() diff --git a/benchmarks/sompy_worker.py b/benchmarks/sompy_worker.py new file mode 100644 index 0000000..4253fab --- /dev/null +++ b/benchmarks/sompy_worker.py @@ -0,0 +1,124 @@ +"""Train one SOMPY map and report how long it took. Runs inside ``.venv-sompy``, not this project. + +SOMPY cannot be imported alongside python-som. It uses ``np.Inf`` at four sites, three of them as +default arguments evaluated when the class body executes, and ``np.Inf`` was removed in NumPy 2.0, +so ``import sompy`` raises ``AttributeError`` on any NumPy this project supports. Two further +reasons to keep it at arm's length even if that were fixed: importing it reconfigures the **root +logger** for the whole process through ``logging.config.dictConfig``, and it pulls in matplotlib and +scikit-image by way of ``sompy.visualization``. + +So the comparison runs across a process boundary. ``bench_vs_sompy.py`` writes a case into a +directory, launches this file with the interpreter from ``.venv-sompy``, and reads the result back. +This half imports **nothing** from python_som, and must not: it would fail on this environment's +NumPy 1.26. + +Not meant to be run by hand. See ``bench_vs_sompy.py`` for the setup command and the protocol. + +Protocol, all inside the directory given as the only argument: + +============== ========================================================================== +``spec.json`` in: ``shape``, ``n_iteration``, ``radius_start``, ``radius_end`` +``input.npz`` in: ``data`` of shape ``(n_samples, n_features)``, ``initial`` of ``(x, y, f)`` +``output.npz`` out: ``weights`` of ``(x, y, f)``, ``seconds`` for the training call alone, + and ``numpy_version``, which the caller cannot read for itself +============== ========================================================================== +""" + +from __future__ import annotations + +import json +import logging +import sys +import time +from pathlib import Path + +import numpy as np +from sompy import SOMFactory + + +def train( + data: np.ndarray, + initial: np.ndarray, + shape: tuple[int, int], + n_iteration: int, + radius_start: float, + radius_end: float, +) -> tuple[np.ndarray, float]: + """Train one SOMPY map from injected models and time the training call alone. + + Three things here are deliberate and each is a control the comparison depends on. + + ``normalization="None"`` selects SOMPY's ``NoNormalizer``, a real pass-through class. At its + ``'var'`` default SOMPY would z-score the data and its models would end up in a different space + from the ones it is being compared against. + + ``_batchtrain`` rather than ``train``, because ``train`` re-runs initialization unconditionally + and would throw the injected codebook away. It is also the only training SOMPY implements: + ``SOMFactory.build`` accepts ``training='seq'`` and then ignores it. + + The codebook is ``(nnodes, n_features)`` with node ``i`` at row ``i // cols``, column + ``i % cols``, so a C-order reshape converts between that and the ``(x, y, f)`` both other + libraries use. + + :param data: Training dataset. + :param initial: Models to start from, of shape ``(x, y, n_features)``. + :param shape: Grid shape. + :param n_iteration: Number of iterations. + :param radius_start: Neighborhood radius at the first iteration. + :param radius_end: Neighborhood radius at the last iteration. + :return: The trained models as ``(x, y, f)``, and the elapsed seconds. + """ + som = SOMFactory.build( + data, + mapsize=list(shape), + normalization="None", + initialization="random", + neighborhood="gaussian", + lattice="rect", + ) + n_features = data.shape[1] + som.codebook.matrix = initial.reshape(som.codebook.nnodes, n_features).copy() + som.codebook.initialized = True + + started = time.perf_counter() + som._batchtrain(trainlen=n_iteration, radiusin=radius_start, radiusfin=radius_end, njob=1) # noqa: SLF001 + elapsed = time.perf_counter() - started + + trained = np.asarray(som.codebook.matrix).reshape(shape[0], shape[1], n_features) + return trained, elapsed + + +def main() -> None: + """Read the case from the directory named on the command line and write the result back.""" + # Before anything else. SOMPY's package __init__ installs a root handler at DEBUG, and + # `_batchtrain` logs an epoch line plus two `timeit` lines per iteration. Left alone that writes + # thousands of lines to a pipe from inside the timed region, which would be measured as SOMPY + # being slow at I/O. `train()` would set the level itself; `_batchtrain` does not. + logging.getLogger().setLevel(logging.ERROR) + + workdir = Path(sys.argv[1]) + spec = json.loads((workdir / "spec.json").read_text(encoding="utf-8")) + with np.load(workdir / "input.npz") as payload: + data, initial = payload["data"], payload["initial"] + + weights, seconds = train( + data, + initial, + (int(spec["shape"][0]), int(spec["shape"][1])), + int(spec["n_iteration"]), + float(spec["radius_start"]), + float(spec["radius_end"]), + ) + # The NumPy version travels with the result. The caller runs a different one by construction, + # so it cannot report this environment's from its own imports, and a comparison that names the + # wrong NumPy beside its numbers is worse than one that names none. + np.savez( + workdir / "output.npz", + weights=weights, + seconds=np.array(seconds), + numpy_version=np.array(np.__version__), + ) + + +if __name__ == "__main__": + main() From 133c14b6473fc2f147227cd27e279ba64ab02653 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Moreira=20Souza?= Date: Thu, 30 Jul 2026 22:45:03 -0300 Subject: [PATCH 4/5] bench: track performance across commits with asv Both benchmark scripts this repository already had compare two implementations at one commit. Neither notices a regression introduced between commits, which is the failure mode a library actually ships. airspeed velocity is what numpy, scipy, pandas and scikit-learn all use for that, and it answers a different question from benchmarks/: this measures python-som against its own git history. Sixteen benchmarks over the parts of this package that have changed before and may change again: both training loops, the neighborhood kernel introduced in 0.4.0 alongside the per-node path it replaced, the matching functions that are the other half of a batch iteration, the SVD that replaced scikit-learn's PCA, and the artifact round trip. Two are peakmem rather than time: the kernel's size is the entire justification for the kernel, since the alternative reached 800 MB on a 100x100 map, and a regression there would otherwise be silent. The directory is named after scikit-learn's so it cannot be confused with benchmarks/, which holds hand-run scripts rather than an asv suite. Nothing gates on a timing. ci.yml gains a step that executes every benchmark once and records nothing, which catches benchmark code rotting unnoticed while leaving the numbers out of the merge decision; the same step checks that the comparison scripts still import and that bench_vs_sompy.py degrades cleanly when its interpreter is absent. Real measurement is workflow_dispatch only, because a shared runner cannot measure anything: an earlier bench_vs_minisom.py swung a comparison from 2.56x to 1.01x between two runs on an idle laptop. Two things about asv's CLI that are not obvious and cost a CI round trip to find out. It has no `check` subcommand, so executing every benchmark without recording is `run --quick --dry-run`. And `asv machine --yes` has to run first: without a machine profile asv exits 1 before running anything. .asv/ is gitignored, for two reasons. Results are per-machine, so a database from one laptop tells nobody anything about another. More importantly, `asv machine` records the host's name, CPU model, core count, RAM and kernel version into .asv/results//machine.json, which is a hardware fingerprint of whoever ran it. The manual workflow uploads only the comparison text for the same reason. asv==0.6.6 goes in a new `bench` extra rather than in `dev`, for the reason `analysis` is separate: it pulls about ten transitive packages and no gating job needs them. --- .github/workflows/benchmarks.yml | 78 ++++++++++++ .github/workflows/ci.yml | 44 +++++++ .gitignore | 8 ++ asv_benchmarks/asv.conf.json | 41 ++++++ asv_benchmarks/benchmarks/__init__.py | 9 ++ asv_benchmarks/benchmarks/artifacts.py | 66 ++++++++++ asv_benchmarks/benchmarks/common.py | 62 +++++++++ asv_benchmarks/benchmarks/initialization.py | 64 ++++++++++ asv_benchmarks/benchmarks/matching.py | 66 ++++++++++ asv_benchmarks/benchmarks/neighborhood.py | 132 ++++++++++++++++++++ asv_benchmarks/benchmarks/training.py | 122 ++++++++++++++++++ pyproject.toml | 15 +++ uv.lock | 127 ++++++++++++++++++- 13 files changed, 833 insertions(+), 1 deletion(-) create mode 100644 .github/workflows/benchmarks.yml create mode 100644 asv_benchmarks/asv.conf.json create mode 100644 asv_benchmarks/benchmarks/__init__.py create mode 100644 asv_benchmarks/benchmarks/artifacts.py create mode 100644 asv_benchmarks/benchmarks/common.py create mode 100644 asv_benchmarks/benchmarks/initialization.py create mode 100644 asv_benchmarks/benchmarks/matching.py create mode 100644 asv_benchmarks/benchmarks/neighborhood.py create mode 100644 asv_benchmarks/benchmarks/training.py diff --git a/.github/workflows/benchmarks.yml b/.github/workflows/benchmarks.yml new file mode 100644 index 0000000..8c39106 --- /dev/null +++ b/.github/workflows/benchmarks.yml @@ -0,0 +1,78 @@ +name: benchmarks + +# Manual only, and that is the point rather than an omission. +# +# asv measures a project against its own git history, and its numbers mean something only as a +# trend on one consistent machine. GitHub's shared runners are not that: an earlier version of +# benchmarks/bench_vs_minisom.py swung a comparison from 2.56x to 1.01x between two runs on an idle +# laptop, and a shared runner is considerably worse than an idle laptop. So nothing here runs on +# push, nothing gates a merge, and the results are an artifact to read rather than a check to pass. +# +# What DOES run on every push is the smoke step in ci.yml, which executes each benchmark once and +# discards the timings. That catches the failure this file cannot: benchmark code rotting unnoticed +# because nobody runs it. +# +# For numbers worth quoting, run it locally on a quiet machine: +# +# cd asv_benchmarks +# uv run --extra bench asv machine --yes +# uv run --extra bench asv continuous master HEAD +on: + workflow_dispatch: + inputs: + base: + description: "Commit, tag or branch to compare against" + required: true + default: master + head: + description: "Commit, tag or branch to measure" + required: true + default: HEAD + factor: + description: "Ratio a change must exceed to be reported" + required: false + default: "1.1" + +# Actions are pinned by commit SHA rather than by tag: a tag can be moved to point at new code, +# a SHA cannot. Dependabot keeps them current. +permissions: + contents: read + +jobs: + asv: + name: asv continuous ${{ inputs.base }}..${{ inputs.head }} + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@3d3c42e5aac5ba805825da76410c181273ba90b1 # v7.0.1 + with: + # asv builds and installs each revision it is asked about, so it needs the history that + # a default shallow clone does not fetch. + fetch-depth: 0 + - uses: astral-sh/setup-uv@c771a70e6277c0a99b617c7a806ffedaca235ff9 # v9.0.0 + with: + enable-cache: true + python-version: "3.12" + - run: uv sync --extra bench + + - name: register the runner + working-directory: asv_benchmarks + run: uv run --extra bench asv machine --yes + + # `|| true` because `asv continuous` exits non-zero when it finds a regression, and a + # regression is the result this job exists to surface rather than a reason to fail it. The + # numbers are in the artifact either way; a human reads them and decides. + - name: compare the two revisions + working-directory: asv_benchmarks + run: | + uv run --extra bench asv continuous \ + --factor "${{ inputs.factor }}" \ + "${{ inputs.base }}" "${{ inputs.head }}" | tee ../asv-comparison.txt || true + + # The comparison text only. `.asv/results/` also contains machine.json, which records the + # host's name, CPU model, core count, RAM and kernel version. That is a hardware fingerprint, + # it is useless to anyone reading the comparison, and a downloadable artifact is a poor place + # to keep it even when the host is an ephemeral runner. + - uses: actions/upload-artifact@043fb46d1a93c77aae656e7c1c64a875d1fc6a0a # v7.0.1 + with: + name: asv-comparison + path: asv-comparison.txt diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index e712728..56e0e7f 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -60,6 +60,50 @@ jobs: - name: pytest run: uv run pytest --cov --cov-report=term-missing --cov-report=xml + benchmarks: + name: benchmarks still run + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@3d3c42e5aac5ba805825da76410c181273ba90b1 # v7.0.1 + - uses: astral-sh/setup-uv@c771a70e6277c0a99b617c7a806ffedaca235ff9 # v9.0.0 + with: + enable-cache: true + python-version: "3.12" + - run: uv sync --extra dev --extra bench + + # Timings are deliberately NOT asserted on. A shared runner cannot measure anything: an + # earlier version of benchmarks/bench_vs_minisom.py swung a comparison from 2.56x to 1.01x + # between two runs on an *idle* machine. What this job protects is that the benchmark code + # still imports and executes, which is the thing that rots silently while nobody runs it. + # + # asv has no `check` subcommand, so the way to execute every benchmark once and record + # nothing is `run --quick --dry-run`. `--python=same` reuses this environment instead of + # building a wheel per revision. `asv machine --yes` is required first: without a machine + # profile asv exits 1 before running anything. + - name: every asv benchmark executes once + working-directory: asv_benchmarks + run: | + uv run --extra bench asv machine --yes + uv run --extra bench asv run --python=same --quick --dry-run + + # The cross-library scripts are not run here. bench_vs_minisom.py trains real maps and takes + # minutes, and bench_vs_sompy.py needs an interpreter this environment cannot provide. Their + # correctness claims are covered by tests/test_minisom_agreement.py, which does run in CI. + # PYTHONPATH rather than sys.path.insert with a relative path, which silently resolves + # against whatever directory the step happens to start in. + - name: the comparison harnesses still import + env: + PYTHONPATH: benchmarks + run: | + uv run python -c " + import bench_update, bench_batch, bench_vs_minisom, bench_vs_sompy + print('benchmark scripts import OK') + " + - name: bench_vs_sompy degrades cleanly without its environment + run: | + test ! -d .venv-sompy + uv run python benchmarks/bench_vs_sompy.py | grep -q "uv venv .venv-sompy" + build: name: build and verify the distribution runs-on: ubuntu-latest diff --git a/.gitignore b/.gitignore index 178b3f2..ab0d8d2 100644 --- a/.gitignore +++ b/.gitignore @@ -113,6 +113,14 @@ venv.bak/ # because it uses np.Inf and so cannot be imported on NumPy 2, so it gets an interpreter of its own. .venv-sompy/ +# Everything airspeed velocity generates: environments, results and the rendered site. +# +# Two reasons, and the second is the stronger one. Results are per-machine, so a database from one +# laptop tells nobody anything about another. And `asv machine` records the host's name, CPU model, +# core count, RAM and kernel version into .asv/results//machine.json, which is a hardware +# fingerprint of whoever ran it and has no business in a public repository. +.asv/ + # Spyder project settings .spyderproject .spyproject diff --git a/asv_benchmarks/asv.conf.json b/asv_benchmarks/asv.conf.json new file mode 100644 index 0000000..7266837 --- /dev/null +++ b/asv_benchmarks/asv.conf.json @@ -0,0 +1,41 @@ +{ + // airspeed velocity configuration. This tracks python-som against its OWN git history; the + // cross-library comparisons live in benchmarks/ and are a different question with a different + // tool. numpy, scipy, pandas and scikit-learn all use asv for exactly this. + // + // Modelled on numpy's, which is the closer fit of the four since this is a pure-Python wheel. + // The directory is named after scikit-learn's so it cannot be confused with benchmarks/, which + // holds hand-run scripts rather than an asv suite. + // + // uv run --extra bench asv run # benchmark the working tree + // uv run --extra bench asv continuous master HEAD # compare two commits + // uv run --extra bench asv publish && asv preview # render the results + // + // Results are per-machine and are NOT committed: a database from one laptop says nothing about + // anyone else's, and .asv/ is gitignored for that reason. + "version": 1, + "project": "python-som", + "project_url": "https://andremsouza.github.io/python-som/", + "repo": "..", + "branches": ["master"], + "dvcs": "git", + "environment_type": "virtualenv", + "show_commit_url": "https://github.com/andremsouza/python-som/commit/", + "build_command": [ + "python -m build --wheel -o {build_cache_dir} {build_dir}" + ], + "install_command": ["python -mpip install {wheel_file}"], + "uninstall_command": ["return-code=any python -mpip uninstall -y {project}"], + "benchmark_dir": "benchmarks", + // Everything asv generates lives under one gitignored directory rather than four. + "env_dir": "../.asv/env", + "results_dir": "../.asv/results", + "html_dir": "../.asv/html", + // NumPy is the only runtime dependency, so an empty list means "whatever resolves", which is + // what a library wanting to notice an upstream regression should track. + "matrix": { + "req": { + "numpy": [] + } + } +} diff --git a/asv_benchmarks/benchmarks/__init__.py b/asv_benchmarks/benchmarks/__init__.py new file mode 100644 index 0000000..598b593 --- /dev/null +++ b/asv_benchmarks/benchmarks/__init__.py @@ -0,0 +1,9 @@ +"""Benchmarks tracking python-som against its own git history. + +Run with ``uv run --extra bench asv run`` from ``asv_benchmarks/``. These measure this package +alone; the cross-library comparisons are hand-run scripts in ``benchmarks/`` and answer a different +question with a different tool. + +Nothing here gates anything. asv numbers are only meaningful as a trend on one consistent machine, +so CI runs each benchmark once to prove it still executes and discards the timings. +""" diff --git a/asv_benchmarks/benchmarks/artifacts.py b/asv_benchmarks/benchmarks/artifacts.py new file mode 100644 index 0000000..1d98a99 --- /dev/null +++ b/asv_benchmarks/benchmarks/artifacts.py @@ -0,0 +1,66 @@ +"""Saving and loading a trained map. + +Not a hot path, and tracked anyway: this is the boundary where a map becomes a file someone keeps, +and the format deliberately refuses pickle. A change that made loading slow enough to discourage +saving would quietly push people back toward ``pickle.dump``, which is arbitrary code execution on +load. Cheap to measure, so it is measured. +""" + +from __future__ import annotations + +import shutil +import tempfile +from pathlib import Path + +import python_som + +from .common import FEATURES, SHAPES, som + + +class Artifacts: + """One save and one load per timed run, on a real temporary file.""" + + params = (SHAPES, FEATURES) + param_names = ("shape", "n_features") + + som: python_som.SOM + directory: Path + path: Path + + def setup(self, shape: tuple[int, int], n_features: int) -> None: + """Build a map and a temporary directory, and write one file to load from. + + :param shape: Grid shape. + :param n_features: Number of features. + """ + self.som = som(shape, n_features) + self.directory = Path(tempfile.mkdtemp(prefix="som-asv-")) + self.path = self.directory / "map.npz" + self.som.save_npz(self.path) + + def teardown(self, shape: tuple[int, int], n_features: int) -> None: + """Remove the temporary directory. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + shutil.rmtree(self.directory, ignore_errors=True) + + def time_save(self, shape: tuple[int, int], n_features: int) -> None: + """Time writing a map, models and metadata together. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + self.som.save_npz(self.directory / "written.npz") + + def time_load(self, shape: tuple[int, int], n_features: int) -> None: + """Time reading one back, including the metadata validation. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + python_som.SOM.load_npz(self.path) diff --git a/asv_benchmarks/benchmarks/common.py b/asv_benchmarks/benchmarks/common.py new file mode 100644 index 0000000..b9ad991 --- /dev/null +++ b/asv_benchmarks/benchmarks/common.py @@ -0,0 +1,62 @@ +"""Shared fixtures and parameter sets for the asv suite. + +Kept in one place so that every benchmark measures the same map on the same data, and so a change +to the case sizes moves all of them together rather than some of them. +""" + +from __future__ import annotations + +import numpy as np + +import python_som + +#: Fixed, so a difference between two commits is a difference in the code. asv reruns the same +#: benchmark against many revisions, and a fresh draw each time would drown the signal it looks for. +SEED = 20260730 + +#: Grid shapes covered by the parameterised benchmarks. Small enough that a full asv sweep over a +#: commit range finishes, large enough that per-node costs are visible. +SHAPES = [(10, 10), (30, 30), (60, 60)] + +#: Feature counts. Varying this shifts how much of the work is the contraction rather than the +#: neighborhood, which is where this package's kernel does or does not pay off. +FEATURES = [4, 16] + +#: Samples in the benchmark dataset. +SAMPLES = 300 + +#: Neighborhood radius. Above the 0.5 default floor, so the floor never engages. +RADIUS = 3.0 + + +def data(n_features: int, n_samples: int = SAMPLES) -> np.ndarray: + """Build the benchmark dataset. + + :param n_features: Number of features. + :param n_samples: Number of samples. + :return: The dataset. + """ + return np.random.default_rng(SEED).normal(size=(n_samples, n_features)) + + +def som(shape: tuple[int, int], n_features: int, neighborhood: str = "gaussian") -> python_som.SOM: + """Build a map with random models, ready to train. + + Random rather than linear initialization: linear would fit a PCA on every setup call, charging + every benchmark for work that ``initialization.py`` measures on its own. + + :param shape: Grid shape. + :param n_features: Number of features. + :param neighborhood: Neighborhood function name. + :return: The map. + """ + built = python_som.SOM( + x=shape[0], + y=shape[1], + input_len=n_features, + neighborhood_radius=RADIUS, + neighborhood_function=neighborhood, + random_seed=SEED, + ) + built.weight_initialization(mode="random") + return built diff --git a/asv_benchmarks/benchmarks/initialization.py b/asv_benchmarks/benchmarks/initialization.py new file mode 100644 index 0000000..b95424b --- /dev/null +++ b/asv_benchmarks/benchmarks/initialization.py @@ -0,0 +1,64 @@ +"""Weight initialization, above all the SVD that replaced scikit-learn's PCA in 0.4.0. + +Dropping scikit-learn saved 264 MB of required install, and the twenty lines of ``np.linalg.svd`` +that replaced it are also more accurate: the solver scikit-learn picked by default was wrong by 5.8% +on data far from the origin. Both claims are checked for correctness by +``tests/test_linalg_matches_sklearn.py``; this is where the cost side stays visible. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING + +import python_som + +from .common import FEATURES, SEED, SHAPES, data + +if TYPE_CHECKING: # pragma: no cover + import numpy as np + + +class Initialize: + """The three ways a map's models can be seeded.""" + + params = (SHAPES, FEATURES) + param_names = ("shape", "n_features") + + som: python_som.SOM + data: np.ndarray + + def setup(self, shape: tuple[int, int], n_features: int) -> None: + """Build an uninitialized map and the dataset outside the timed region. + + :param shape: Grid shape. + :param n_features: Number of features. + """ + self.som = python_som.SOM(x=shape[0], y=shape[1], input_len=n_features, random_seed=SEED) + self.data = data(n_features) + + def time_linear(self, shape: tuple[int, int], n_features: int) -> None: + """Time the PCA-based initialization, which is the SVD. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + self.som.weight_initialization(mode="linear", data=self.data) + + def time_sample(self, shape: tuple[int, int], n_features: int) -> None: + """Time seeding from drawn samples. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + self.som.weight_initialization(mode="sample", data=self.data) + + def time_random(self, shape: tuple[int, int], n_features: int) -> None: + """Time the cheapest initializer, as a floor to read the other two against. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + self.som.weight_initialization(mode="random") diff --git a/asv_benchmarks/benchmarks/matching.py b/asv_benchmarks/benchmarks/matching.py new file mode 100644 index 0000000..f22413c --- /dev/null +++ b/asv_benchmarks/benchmarks/matching.py @@ -0,0 +1,66 @@ +"""Matching inputs to models: the other half of batch training, and all of scoring. + +``accumulate`` is roughly 40% of a batch iteration and loops over samples in Python, so it is the +half of the cost that scales with the dataset rather than with the grid. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING + +from python_som._core._distance import euclidean_distance +from python_som._core._match import accumulate, quantization, winner + +from .common import FEATURES, SHAPES, data, som + +if TYPE_CHECKING: # pragma: no cover + import numpy as np + + +class Match: + """Winner search, accumulation and quantization on the same map and dataset.""" + + params = (SHAPES, FEATURES) + param_names = ("shape", "n_features") + + weights: np.ndarray + data: np.ndarray + shape: tuple[int, int] + + def setup(self, shape: tuple[int, int], n_features: int) -> None: + """Build the models and dataset outside the timed region. + + :param shape: Grid shape. + :param n_features: Number of features. + """ + self.shape = shape + self.weights = som(shape, n_features).get_weights() + self.data = data(n_features) + + def time_accumulate(self, shape: tuple[int, int], n_features: int) -> None: + """Time one batch iteration's worth of per-node sums and counts. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + accumulate(self.data, self.weights, self.shape, euclidean_distance) + + def time_winner_for_every_sample(self, shape: tuple[int, int], n_features: int) -> None: + """Time the best-matching-unit search alone, without the accumulation around it. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + for sample in self.data: + winner(sample, self.weights, euclidean_distance) + + def time_quantization(self, shape: tuple[int, int], n_features: int) -> None: + """Time what ``quantization_error`` costs, which every ``train`` call pays once. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + quantization(self.data, self.weights, euclidean_distance) diff --git a/asv_benchmarks/benchmarks/neighborhood.py b/asv_benchmarks/benchmarks/neighborhood.py new file mode 100644 index 0000000..42b9e4b --- /dev/null +++ b/asv_benchmarks/benchmarks/neighborhood.py @@ -0,0 +1,132 @@ +"""The neighborhood kernel, which is the optimization 0.4.0 rests on. + +Batch training needs ``h_ji`` for every pair of nodes. Because a neighborhood depends only on the +offset between two nodes, one kernel over every offset serves the whole grid and each node's +neighborhood is a slice of it. Evaluating per node instead was 42% of batch training on a 40x40 map. + +Both paths are still benchmarked, and the per-node one is the point: it is what the kernel +replaced, so keeping it measured is what makes the claim checkable rather than historical. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING + +from python_som._core._neighborhood import kernel_view, resolve, resolve_kernel + +from .common import RADIUS, SHAPES + +if TYPE_CHECKING: # pragma: no cover + from collections.abc import Callable + + import numpy as np + import numpy.typing as npt + +#: Both neighborhoods batch training accepts, plus the signed one only stepwise can use. +NEIGHBORHOODS = ["gaussian", "bubble", "mexican_hat"] + +#: No wrapping, wrapping on one axis, wrapping on both. The cyclic fold is the part of the offset +#: machinery most likely to be got wrong, and the part whose cost is least obvious. +CYCLIC = [(False, False), (True, False), (True, True)] + + +class Kernel: + """Building one kernel per iteration.""" + + params = (SHAPES, NEIGHBORHOODS, CYCLIC) + param_names = ("shape", "neighborhood", "cyclic") + + build: Callable[..., npt.NDArray[np.floating]] + + def setup(self, shape: tuple[int, int], neighborhood: str, cyclic: tuple[bool, bool]) -> None: + """Resolve the kernel builder outside the timed region. + + :param shape: Grid shape. + :param neighborhood: Neighborhood function name. + :param cyclic: Whether each axis wraps. + """ + del shape, cyclic + self.build = resolve_kernel(neighborhood) + + def time_build( + self, shape: tuple[int, int], neighborhood: str, cyclic: tuple[bool, bool] + ) -> None: + """Time building the kernel once. + + :param shape: Grid shape. + :param neighborhood: Unused. + :param cyclic: Whether each axis wraps. + """ + del neighborhood + self.build(shape, RADIUS, cyclic) + + def peakmem_build( + self, shape: tuple[int, int], neighborhood: str, cyclic: tuple[bool, bool] + ) -> None: + """Track the kernel's size, which is the justification for the whole approach. + + A ``(2X-1, 2Y-1)`` kernel is 198 KB at 80x80 against the 800 MB a full ``(x, y, x, y)`` + tensor would need. A regression here would otherwise be silent. + + :param shape: Grid shape. + :param neighborhood: Unused. + :param cyclic: Whether each axis wraps. + """ + del neighborhood + self.build(shape, RADIUS, cyclic) + + +class Slice: + """Taking one node's neighborhood out of a built kernel. + + Must stay a view rather than a copy. Copying ``(X, Y)`` floats per node would give back most of + what the kernel wins, and this is where that would show up as a trend. + """ + + params = (SHAPES,) + param_names = ("shape",) + + kernel: npt.NDArray[np.floating] + nodes: list[tuple[int, int]] + + def setup(self, shape: tuple[int, int]) -> None: + """Build the kernel and the node list outside the timed region. + + :param shape: Grid shape. + """ + self.kernel = resolve_kernel("gaussian")(shape, RADIUS, (False, False)) + self.nodes = [(x, y) for x in range(shape[0]) for y in range(shape[1])] + + def time_slice_every_node(self, shape: tuple[int, int]) -> None: + """Time slicing the kernel once per node, which is one batch iteration's worth. + + :param shape: Grid shape. + """ + for node in self.nodes: + kernel_view(self.kernel, shape, node) + + +class PerNode: + """Evaluating the neighborhood once per node, which the kernel replaced.""" + + params = (SHAPES,) + param_names = ("shape",) + + evaluate: Callable[..., npt.NDArray[np.floating]] + nodes: list[tuple[int, int]] + + def setup(self, shape: tuple[int, int]) -> None: + """Resolve the per-node function and build the node list outside the timed region. + + :param shape: Grid shape. + """ + self.evaluate = resolve("gaussian") + self.nodes = [(x, y) for x in range(shape[0]) for y in range(shape[1])] + + def time_evaluate_every_node(self, shape: tuple[int, int]) -> None: + """Time the path the kernel replaced, on the same work. + + :param shape: Grid shape. + """ + for node in self.nodes: + self.evaluate(shape, node, RADIUS, (False, False)) diff --git a/asv_benchmarks/benchmarks/training.py b/asv_benchmarks/benchmarks/training.py new file mode 100644 index 0000000..956615d --- /dev/null +++ b/asv_benchmarks/benchmarks/training.py @@ -0,0 +1,122 @@ +"""The training loops, which are what a performance change in this package almost always touches.""" + +from __future__ import annotations + +from typing import TYPE_CHECKING, Any + +from .common import FEATURES, SHAPES, data, som + +if TYPE_CHECKING: # pragma: no cover + import numpy as np + + import python_som + +#: Iterations per timed run. Batch is the more expensive mode per iteration, so it gets fewer. +BATCH_ITERATIONS = 10 +STEPWISE_ITERATIONS = 200 + + +class Batch: + """Kohonen Eq. (8), the mode this package treats as primary.""" + + params = (SHAPES, FEATURES) + param_names = ("shape", "n_features") + + som: python_som.SOM + data: np.ndarray + + def setup(self, shape: tuple[int, int], n_features: int) -> None: + """Build the map and dataset outside the timed region. + + :param shape: Grid shape. + :param n_features: Number of features. + """ + self.som = som(shape, n_features) + self.data = data(n_features) + + def time_train(self, shape: tuple[int, int], n_features: int) -> None: + """Time batch training. + + The private loop rather than ``train``, which would also score the whole dataset to fill in + its ``TrainingReport``. That pass is real work a user pays for, but it is not training, and + folding it in here would blur a change in the loop with a change in the metric. + + :param shape: Unused; asv passes the parameters back. + :param n_features: Unused. + """ + del shape, n_features + self.som._train_batch(self.data, BATCH_ITERATIONS, verbose=False) # noqa: SLF001 + + def peakmem_train(self, shape: tuple[int, int], n_features: int) -> None: + """Track the memory batch training holds at once. + + Worth tracking rather than assuming: the neighborhood kernel this package introduced in + 0.4.0 is a deliberate memory decision, and the alternative it replaced was a tensor that + reached 800 MB on a 100x100 map. + + :param shape: Unused. + :param n_features: Unused. + """ + del shape, n_features + self.som._train_batch(self.data, BATCH_ITERATIONS, verbose=False) # noqa: SLF001 + + +class Stepwise: + """Kohonen Eq. (3), one sample at a time.""" + + params = (SHAPES, FEATURES, ["random", "sequential"]) + param_names = ("shape", "n_features", "mode") + + som: python_som.SOM + data: np.ndarray + + def setup(self, shape: tuple[int, int], n_features: int, mode: str) -> None: + """Build the map and dataset outside the timed region. + + :param shape: Grid shape. + :param n_features: Number of features. + :param mode: Training mode. + """ + del mode + self.som = som(shape, n_features) + self.data = data(n_features) + + def time_train(self, shape: tuple[int, int], n_features: int, mode: str) -> None: + """Time stepwise training in either sample order. + + :param shape: Unused. + :param n_features: Unused. + :param mode: Either ``'random'`` or ``'sequential'``. + """ + del shape, n_features + self.som._train_stepwise( # noqa: SLF001 + self.data, STEPWISE_ITERATIONS, mode=mode, verbose=False + ) + + +class MexicanHat: + """The signed neighborhood, which batch rejects and so only stepwise can exercise.""" + + params = (SHAPES,) + param_names = ("shape",) + + som: python_som.SOM + data: Any + + def setup(self, shape: tuple[int, int]) -> None: + """Build the map and dataset outside the timed region. + + :param shape: Grid shape. + """ + self.som = som(shape, FEATURES[0], neighborhood="mexican_hat") + self.data = data(FEATURES[0]) + + def time_train(self, shape: tuple[int, int]) -> None: + """Time stepwise training with the mexican hat. + + :param shape: Unused. + """ + del shape + self.som._train_stepwise( # noqa: SLF001 + self.data, STEPWISE_ITERATIONS, mode="random", verbose=False + ) diff --git a/pyproject.toml b/pyproject.toml index 9937bed..8bf3e98 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -102,6 +102,21 @@ docs = [ # nothing in `src/`, and this is a pure-numpy library with no network, auth, crypto, SQL or # templating, which is most of what semgrep's rulesets look for. Low coverage gained for three CVEs # that cannot be patched. Revisit if semgrep relaxes the pin. +# Benchmarking, kept out of `dev` for the same reason `analysis` is: `asv` pulls about ten +# transitive packages (asv-runner, json5, build, tabulate, virtualenv, packaging, +# importlib-metadata, pyyaml, pympler) and no gating CI job runs it. `minisom` is repeated here so +# that `--extra bench` alone is enough to run every script in benchmarks/. +# +# SOMPY is deliberately absent and cannot be added. It uses `np.Inf`, removed in NumPy 2.0, at +# class-definition time, so it cannot be imported in any environment this package supports. +# `benchmarks/bench_vs_sompy.py` drives it through a separate interpreter instead; that environment +# is built by hand and the command is in the script's docstring. +bench = [ + # airspeed velocity, from the airspeed-velocity organisation: 22 releases since 2015-05-01, + # and the tool numpy, scipy, pandas and scikit-learn all use to track their own performance. + "asv==0.6.6", + "minisom==2.3.6", +] analysis = [ # PyCQA, first released 2015. Security-oriented AST checks. "bandit==1.9.4", diff --git a/uv.lock b/uv.lock index 2da1ee0..cf37b41 100644 --- a/uv.lock +++ b/uv.lock @@ -70,6 +70,45 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/b0/cf/1c5f42b110e57bc5502eb80dbc3b03d256926062519224835ef08134f1f9/astroid-4.0.4-py3-none-any.whl", hash = "sha256:52f39653876c7dec3e3afd4c2696920e05c83832b9737afc21928f2d2eb7a753", size = 276445, upload-time = "2026-02-07T23:35:05.344Z" }, ] +[[package]] +name = "asv" +version = "0.6.6" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "asv-runner" }, + { name = "build" }, + { name = "colorama", marker = "sys_platform == 'win32'" }, + { name = "importlib-metadata" }, + { name = "json5" }, + { name = "packaging" }, + { name = "pympler", marker = "platform_python_implementation != 'PyPy'" }, + { name = "pyyaml", marker = "platform_python_implementation != 'PyPy'" }, + { name = "tabulate" }, + { name = "tomli", marker = "python_full_version < '3.11'" }, + { name = "virtualenv" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/a4/0d/17a15f8d941c6846d06708d33ae03c4154844fd7c9610f12be6837241b05/asv-0.6.6.tar.gz", hash = "sha256:82e47105db8f56d9b1e54763dd01a1709d2722ad2f727b21521628c3e110bdc6", size = 403111, upload-time = "2026-06-27T03:12:47.749Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/5f/dc/364295d8488b8b325efeb57daf6c121d32749a489cacca7446251525d7e4/asv-0.6.6-cp314-cp314t-macosx_10_15_x86_64.whl", hash = "sha256:7f66ceff065fa02c342a00ccf9832ec34dca3835493173e9cf199851f6686c2b", size = 187438, upload-time = "2026-06-27T03:12:26.272Z" }, + { url = "https://files.pythonhosted.org/packages/35/65/e8110d7e833a3ff44fbbf00d9267025cf6c31d2d28968897a1dffd149a6c/asv-0.6.6-cp314-cp314t-macosx_11_0_arm64.whl", hash = "sha256:68bdabaf4c4441c460dfe2b9c1722a7f24f0c5cc2f284a751a3fbee75882c87a", size = 187995, upload-time = "2026-06-27T03:12:28.08Z" }, + { url = "https://files.pythonhosted.org/packages/c9/e8/35a60d4c0dc90b319ed80e61cba49668c1c92a5ab289ca2062b3764ae8e3/asv-0.6.6-cp314-cp314t-manylinux2014_x86_64.manylinux_2_17_x86_64.whl", hash = "sha256:dfdc4a6295c8539be8c11136d7aaabc6e4293efbc9f635f0263d02205bdd53a2", size = 293614, upload-time = "2026-06-27T03:13:15.432Z" }, + { url = "https://files.pythonhosted.org/packages/b3/1c/98cc990d856e84cd5f1212e978fcc4d27e48fbc97707a222d295492812d5/asv-0.6.6-cp314-cp314t-musllinux_1_2_x86_64.whl", hash = "sha256:e2d47388069730ded8c0955fdeb8369d810641faa44c3687864ef8941782f2d1", size = 1315341, upload-time = "2026-06-27T03:13:21.731Z" }, + { url = "https://files.pythonhosted.org/packages/ec/0c/a6e0e9b9b7c84eab45001a1cd0632634405fe7c8c7332f315b727cc87cf3/asv-0.6.6-cp314-cp314t-win_amd64.whl", hash = "sha256:acfaf32d34301bd1b7386d533f4005b5f94b7cf0c0509042ee942f02ff4de3d5", size = 189115, upload-time = "2026-06-27T03:13:41.891Z" }, + { url = "https://files.pythonhosted.org/packages/59/fc/83cf0d5ae0052d7ebc8c5e08f6ed2c2bd15c5d895b82652ec09cf21fa317/asv-0.6.6-cp36-abi3-macosx_10_9_x86_64.whl", hash = "sha256:061cd2c370b3427ccf4bdab7c9a6f7ca593b7b74f6f463b825809840b73371a2", size = 186451, upload-time = "2026-06-27T03:12:29.55Z" }, + { url = "https://files.pythonhosted.org/packages/be/34/3ba6fd0fd33de36238a185f6192711f7c169d164e5ad3895bd0baa9414c1/asv-0.6.6-cp36-abi3-macosx_11_0_arm64.whl", hash = "sha256:a4a70ad4a4cd45c7e6d72f1febc56c0b092ccc21fd06132e9806413a336a97d7", size = 187146, upload-time = "2026-06-27T03:12:31.042Z" }, + { url = "https://files.pythonhosted.org/packages/70/be/b35663e3ffe0ee23b4d1be754816b1d676d9a65418220d675adb7502db06/asv-0.6.6-cp36-abi3-manylinux2014_x86_64.manylinux_2_17_x86_64.whl", hash = "sha256:0272503beb40b21fbdeb9149b290275791fd824a3a297ffa5fb7029ace989636", size = 280305, upload-time = "2026-06-27T03:13:28.48Z" }, + { url = "https://files.pythonhosted.org/packages/98/c6/1d96bc66440bfe2156020b9655522dfc3839ff4137d72f09a44a131d0a91/asv-0.6.6-cp36-abi3-musllinux_1_2_x86_64.whl", hash = "sha256:93e6480d87965a60573fe9d48645860f9cd4958a7bfb3b080f43ec08a96e62ee", size = 1303571, upload-time = "2026-06-27T03:13:34.675Z" }, + { url = "https://files.pythonhosted.org/packages/e0/21/7106ffc89c07faad7b56af3d7166cbc6a25c21845695365be1a06ce44b1e/asv-0.6.6-cp36-abi3-win_amd64.whl", hash = "sha256:a18a2bf9441bfe55f34f0f192db178eee6ead219e11752732f4e9230353fa8d4", size = 188189, upload-time = "2026-06-27T03:13:47.955Z" }, +] + +[[package]] +name = "asv-runner" +version = "0.2.5" +source = { registry = "https://pypi.org/simple" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/4d/30/8d4ac65662f20dd8e9c662ed1ca7c1ad11bcf39e40a776b4b837ab7adc86/asv_runner-0.2.5-py3-none-any.whl", hash = "sha256:89b7d68018bae3629022fb2bfb065699f12ed2e5a7022bd6f85110fc713463b1", size = 50577, upload-time = "2026-06-27T01:41:04.712Z" }, +] + [[package]] name = "babel" version = "2.18.0" @@ -126,6 +165,22 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/e5/ca/78d423b324b8d77900030fa59c4aa9054261ef0925631cd2501dd015b7b7/boolean_py-5.0-py3-none-any.whl", hash = "sha256:ef28a70bd43115208441b53a045d1549e2f0ec6e3d08a9d142cbc41c1938e8d9", size = 26577, upload-time = "2025-04-03T10:39:48.449Z" }, ] +[[package]] +name = "build" +version = "1.5.0" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "colorama", marker = "os_name == 'nt'" }, + { name = "importlib-metadata", marker = "python_full_version < '3.10.2'" }, + { name = "packaging" }, + { name = "pyproject-hooks" }, + { name = "tomli", marker = "python_full_version < '3.11'" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/78/e0/df5e171f685f82f37b12e1f208064e24244911079d7b767447d1af7e0d70/build-1.5.0.tar.gz", hash = "sha256:302c22c3ba2a0fd5f3911918651341ebb3896176cbdec15bd421f80b1afc7647", size = 89796, upload-time = "2026-04-30T03:18:25.17Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/0d/fe/6bea5c9162869c5beba5d9c8abbed835ec85bf1ec1fba05a3822325c45f3/build-1.5.0-py3-none-any.whl", hash = "sha256:13f3eecb844759ab66efec90ca17639bbf14dc06cb2fdf37a9010322d9c50a6f", size = 26018, upload-time = "2026-04-30T03:18:23.644Z" }, +] + [[package]] name = "cachecontrol" version = "0.14.4" @@ -994,6 +1049,15 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/7b/91/984aca2ec129e2757d1e4e3c81c3fcda9d0f85b74670a094cc443d9ee949/joblib-1.5.3-py3-none-any.whl", hash = "sha256:5fc3c5039fc5ca8c0276333a188bbd59d6b7ab37fe6632daa76bc7f9ec18e713", size = 309071, upload-time = "2025-12-15T08:41:44.973Z" }, ] +[[package]] +name = "json5" +version = "0.15.0" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/e4/7d/05c46a96a78147ae3bf99c2f4169ce144a70220b8d6fcd56f6ec368b8ce9/json5-0.15.0.tar.gz", hash = "sha256:7424d1f1eb1d56da6e3d70643f53619862b4ce81440bdb8ecfd6f875e5ba4a71", size = 53278, upload-time = "2026-06-19T20:08:27.716Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/eb/be/59527c99478aade6bb33a68d72e6e18dd4e6ff6eacfc7d01bdb15bc76912/json5-0.15.0-py3-none-any.whl", hash = "sha256:56636a30c0e8a4665fe2179c0212f32eae3796dea89ea6f649b9436ecdb39618", size = 36570, upload-time = "2026-06-19T20:08:26.748Z" }, +] + [[package]] name = "keyring" version = "25.7.0" @@ -2497,6 +2561,18 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/d6/54/da572c98c0b77626a91b5d3b89f0231d8bff5125c225420908632f8b342d/pymdown_extensions-11.0.1-py3-none-any.whl", hash = "sha256:db3943a62bab7e03af1364f0c4083e64b91fb097675a4b6cceccfbe9a77e5eb2", size = 269455, upload-time = "2026-07-02T17:59:21.271Z" }, ] +[[package]] +name = "pympler" +version = "1.1" +source = { registry = "https://pypi.org/simple" } +dependencies = [ + { name = "pywin32", marker = "sys_platform == 'win32'" }, +] +sdist = { url = "https://files.pythonhosted.org/packages/dd/37/c384631908029676d8e7213dd956bb686af303a80db7afbc9be36bc49495/pympler-1.1.tar.gz", hash = "sha256:1eaa867cb8992c218430f1708fdaccda53df064144d1c5656b1e6f1ee6000424", size = 179954, upload-time = "2024-06-28T19:56:06.563Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/79/4f/a6a2e2b202d7fd97eadfe90979845b8706676b41cbd3b42ba75adf329d1f/Pympler-1.1-py3-none-any.whl", hash = "sha256:5b223d6027d0619584116a0cbc28e8d2e378f7a79c1e5e024f9ff3b673c58506", size = 165766, upload-time = "2024-06-28T19:56:05.087Z" }, +] + [[package]] name = "pyparsing" version = "3.3.2" @@ -2506,6 +2582,15 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/10/bd/c038d7cc38edc1aa5bf91ab8068b63d4308c66c4c8bb3cbba7dfbc049f9c/pyparsing-3.3.2-py3-none-any.whl", hash = "sha256:850ba148bd908d7e2411587e247a1e4f0327839c40e2e5e6d05a007ecc69911d", size = 122781, upload-time = "2026-01-21T03:57:55.912Z" }, ] +[[package]] +name = "pyproject-hooks" +version = "1.2.0" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/e7/82/28175b2414effca1cdac8dc99f76d660e7a4fb0ceefa4b4ab8f5f6742925/pyproject_hooks-1.2.0.tar.gz", hash = "sha256:1e859bd5c40fae9448642dd871adf459e5e2084186e8d2c2a79a824c970da1f8", size = 19228, upload-time = "2024-09-29T09:24:13.293Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/bd/24/12818598c362d7f300f18e74db45963dbcb85150324092410c8b49405e42/pyproject_hooks-1.2.0-py3-none-any.whl", hash = "sha256:9e5c6bfa8dcc30091c74b0cf803c81fdd29d94f01992a7707bc97babb1141913", size = 10216, upload-time = "2024-09-29T09:24:11.978Z" }, +] + [[package]] name = "pytest" version = "9.1.1" @@ -2579,6 +2664,10 @@ analysis = [ { name = "pip-audit" }, { name = "pylint" }, ] +bench = [ + { name = "asv" }, + { name = "minisom" }, +] cli = [ { name = "tqdm" }, ] @@ -2620,9 +2709,11 @@ sklearn = [ [package.metadata] requires-dist = [ + { name = "asv", marker = "extra == 'bench'", specifier = "==0.6.6" }, { name = "bandit", marker = "extra == 'analysis'", specifier = "==1.9.4" }, { name = "hypothesis", marker = "extra == 'dev'", specifier = "==6.163.0" }, { name = "matplotlib", marker = "extra == 'examples'", specifier = ">=3.8" }, + { name = "minisom", marker = "extra == 'bench'", specifier = "==2.3.6" }, { name = "minisom", marker = "extra == 'dev'", specifier = "==2.3.6" }, { name = "mkdocs-material", marker = "extra == 'docs'", specifier = "==9.7.7" }, { name = "mkdocs-redirects", marker = "extra == 'docs'", specifier = "==1.2.2" }, @@ -2649,7 +2740,7 @@ requires-dist = [ { name = "twine", marker = "extra == 'dev'", specifier = "==7.0.0" }, { name = "types-tqdm", marker = "extra == 'dev'", specifier = "==4.69.0.20260728" }, ] -provides-extras = ["cli", "dev", "docs", "analysis", "sklearn", "examples"] +provides-extras = ["cli", "dev", "docs", "bench", "analysis", "sklearn", "examples"] [[package]] name = "pytz" @@ -2660,6 +2751,31 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/0f/7b/39c34ca613b0b198cb866466651b26b045e2009864c5183c979a3b83f383/pytz-2026.3.post1-py2.py3-none-any.whl", hash = "sha256:dd95840dd199baea12d9cc096a1d452caa6596a1c1e4b5f3dbd1541855d5e815", size = 508283, upload-time = "2026-07-25T15:12:05.782Z" }, ] +[[package]] +name = "pywin32" +version = "312" +source = { registry = "https://pypi.org/simple" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/fe/1b/9cfdeac80ee45bebbbcb31f1b7b99a0d81a1c72de48d837be984e0e88b1d/pywin32-312-cp310-cp310-win32.whl", hash = "sha256:772235332b5d1024c696f11cea1ae4be7930f0a8b894bb43db14e3f435f1ff7e", size = 6361387, upload-time = "2026-06-04T07:49:14.329Z" }, + { url = "https://files.pythonhosted.org/packages/33/b1/7afc96d041d982c27bc2df6f853d43f01fd273e3d39d04be3647ddeb533d/pywin32-312-cp310-cp310-win_amd64.whl", hash = "sha256:5dbc35d2b5320dc07f25fa31269cfb767471002b17de5eb067d03da68c7cb2db", size = 6926780, upload-time = "2026-06-04T07:49:16.881Z" }, + { url = "https://files.pythonhosted.org/packages/ce/3a/4140da9ad54108e517f4a16b2d83da3033e08662144623e1239587cb7db6/pywin32-312-cp310-cp310-win_arm64.whl", hash = "sha256:3020656e34f1cf7faeb7bccd2b84653a607c6ff0c55ada85e6487d61716deabd", size = 4307203, upload-time = "2026-06-04T07:49:18.993Z" }, + { url = "https://files.pythonhosted.org/packages/1f/f5/10a6e845a00fc5e7afd0a988b744f403d4d57162a28d160a093c4d9322f0/pywin32-312-cp311-cp311-win32.whl", hash = "sha256:17948aeadbdb091f0ced6ef0841620794e68327b94ee415571c1203594b7215c", size = 6362659, upload-time = "2026-06-04T07:49:21.349Z" }, + { url = "https://files.pythonhosted.org/packages/35/c4/dcd2d62b5944b6d5db53413a5899016ccd57ffcb7278f3f81655d25d2027/pywin32-312-cp311-cp311-win_amd64.whl", hash = "sha256:d11417d84412f859b722fad0841b3614459ed0047f7542d8362e77884f6b6e8a", size = 6928825, upload-time = "2026-06-04T07:49:23.934Z" }, + { url = "https://files.pythonhosted.org/packages/b7/56/3cbb433fe4501cdba2eb9040f56a4e1a8243faa4186b25295564d1a7a79d/pywin32-312-cp311-cp311-win_arm64.whl", hash = "sha256:b2200a054ca6d6625c4842fc56a4976a4b47f96b73dbe5538c3f813a80359f47", size = 6721875, upload-time = "2026-06-04T07:49:26.416Z" }, + { url = "https://files.pythonhosted.org/packages/83/ff/32aa7d2ed0ab12b323aaa64f9b75e6ad4f8fd09f9ccfc28c79414d46838d/pywin32-312-cp312-cp312-win32.whl", hash = "sha256:dab4f65ac9c4e48400a2a0530c46c3c579cd5905ecd11b80692373915269208b", size = 6371877, upload-time = "2026-06-04T07:49:28.836Z" }, + { url = "https://files.pythonhosted.org/packages/03/d9/77040d3b43df3f3be32ea289433d660d2727f5ba327bc73be835127d9d60/pywin32-312-cp312-cp312-win_amd64.whl", hash = "sha256:b457f6d628a47e8a7346ce22acb7e1a46a4a78b52e1d17e1af56871bd19a93bc", size = 6914841, upload-time = "2026-06-04T07:49:31.85Z" }, + { url = "https://files.pythonhosted.org/packages/e3/cc/7b1ec671775756020a0ee7f4feeaf3c568f0ab86bd3900088cf986937a92/pywin32-312-cp312-cp312-win_arm64.whl", hash = "sha256:6017c58e12f6809fbb0555b75df144c2922a9ffd18e4b9b5afa863b6c1a9d950", size = 6727901, upload-time = "2026-06-04T07:49:34.244Z" }, + { url = "https://files.pythonhosted.org/packages/2d/41/12fbfd7f36ed2146d8bc9de96c2741296bf0d490b98508496cff322e274c/pywin32-312-cp313-cp313-win32.whl", hash = "sha256:7a27df850933d16a8eabfbaeb73d52b273e2da667f80d70b01a89d1f6828d02c", size = 6370184, upload-time = "2026-06-04T07:49:36.253Z" }, + { url = "https://files.pythonhosted.org/packages/ba/db/36a78e3403099d31d9746d13fdcde5accc43c1155f375a34d15983a479a7/pywin32-312-cp313-cp313-win_amd64.whl", hash = "sha256:c53e878d15a1c44788082bfe712a905433473aa38f86375b7cf8b45e3acbaaf9", size = 6914298, upload-time = "2026-06-04T07:49:38.876Z" }, + { url = "https://files.pythonhosted.org/packages/84/37/c1697194092b76de9ed47ca124323f02c57ffc8a45c06f88a3d5acaf01eb/pywin32-312-cp313-cp313-win_arm64.whl", hash = "sha256:59aba5d5940842075343a5ddc6b11f1cdf0d1567fe745290359dfbcc7c2eb831", size = 6727640, upload-time = "2026-06-04T07:49:41.083Z" }, + { url = "https://files.pythonhosted.org/packages/fc/2b/1f3cded5822fd49c02f40544cbb5f58c7cfd6b1694869fd476cb6170ee97/pywin32-312-cp314-cp314-win32.whl", hash = "sha256:a77a90fbb6881238d2ca9c6fd797b25817f3768fe78d214a90137ff055a75f5b", size = 6468928, upload-time = "2026-06-04T07:49:43.188Z" }, + { url = "https://files.pythonhosted.org/packages/21/82/3bf86d2e2808902013132e1ce905a7da0da53790f3836c64bf44d55e24f3/pywin32-312-cp314-cp314-win_amd64.whl", hash = "sha256:a4dd3a848290ef724347b19f301045831d8e802fa4464f491b98b1e0a081432e", size = 7024157, upload-time = "2026-06-04T07:49:45.34Z" }, + { url = "https://files.pythonhosted.org/packages/a4/0e/73f6d6800b4f27655abd9e9f6aaeaefcddb2b946e4674efa2bab184a7f7b/pywin32-312-cp314-cp314-win_arm64.whl", hash = "sha256:9fce94568364e0155e6dfb781ac5d95903be8baf28670632beab1b523f300daa", size = 6839598, upload-time = "2026-06-04T07:49:47.613Z" }, + { url = "https://files.pythonhosted.org/packages/eb/61/caa39686032d2ebdd04ff0ab5cbe163126c0066d98e00c9018646e42393b/pywin32-312-cp315-cp315-win32.whl", hash = "sha256:5c1fbe4a937a73ae9297384a3da38518cbc694c68ad8a809b2e19acd350f03ed", size = 6471159, upload-time = "2026-06-04T07:49:50.035Z" }, + { url = "https://files.pythonhosted.org/packages/0f/cd/7e1de64a4a6f69c04214169657ccab0d93a670ea50e35eb8f489d7378249/pywin32-312-cp315-cp315-win_amd64.whl", hash = "sha256:c2f03a0f73f804a13c2735b99392b0cd426bb4f2c4d0178e5ac966a0f21618d5", size = 7025293, upload-time = "2026-06-04T07:49:54.857Z" }, + { url = "https://files.pythonhosted.org/packages/23/ed/4532e9388e65fa16b46776ef47ad631a64eda1631884488af707666350ed/pywin32-312-cp315-cp315-win_arm64.whl", hash = "sha256:a8597d28f267b39074aef51fa593530082b39cbe5a074226096857b1fed2dfb9", size = 6840337, upload-time = "2026-06-04T07:49:57.531Z" }, +] + [[package]] name = "pywin32-ctypes" version = "0.2.3" @@ -3222,6 +3338,15 @@ wheels = [ { url = "https://files.pythonhosted.org/packages/62/8d/008761f6e1000600e5303db30d05724bdcf3d2d186cbb59fac79b52e39ed/stevedore-5.9.0-py3-none-any.whl", hash = "sha256:e520945d4c257700eddc1eb1d79df04b2ea578eef185e0e3fa5b442fc848d3f7", size = 54463, upload-time = "2026-07-02T11:38:07.43Z" }, ] +[[package]] +name = "tabulate" +version = "0.10.0" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/46/58/8c37dea7bbf769b20d58e7ace7e5edfe65b849442b00ffcdd56be88697c6/tabulate-0.10.0.tar.gz", hash = "sha256:e2cfde8f79420f6deeffdeda9aaec3b6bc5abce947655d17ac662b126e48a60d", size = 91754, upload-time = "2026-03-04T18:55:34.402Z" } +wheels = [ + { url = "https://files.pythonhosted.org/packages/99/55/db07de81b5c630da5cbf5c7df646580ca26dfaefa593667fc6f2fe016d2e/tabulate-0.10.0-py3-none-any.whl", hash = "sha256:f0b0622e567335c8fabaaa659f1b33bcb6ddfe2e496071b743aa113f8774f2d3", size = 39814, upload-time = "2026-03-04T18:55:31.284Z" }, +] + [[package]] name = "threadpoolctl" version = "3.6.0" From 7df0ba0c59218e8d6e26489b405fdc115411d672 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andr=C3=A9=20Moreira=20Souza?= Date: Thu, 30 Jul 2026 22:51:00 -0300 Subject: [PATCH 5/5] docs: publish the comparison with MiniSom and SOMPY The measurements exist in benchmarks/ but nobody choosing between these three libraries reads a benchmark script. This is the page that answers the question, and most of its length goes to why the measurement is harder than it looks rather than to the numbers. Both results are published as measured. Against SOMPY this package is 1.3x to 2.0x faster; against MiniSom it is 1.3x to 1.6x SLOWER at batch training, with the gap growing as the map grows, and the page says why: batch_update walks every node in Python where MiniSom's node-side update is a single vectorised divide, so the ratio tracks nodes against samples. Stepwise is within 10%. The page leads with the differences rather than the timings, because three packages that look interchangeable are not: MiniSom's train_batch is stepwise training, only the gaussian is the same function in all three, and the bubble selects 9 nodes here against 1 in MiniSom at sigma=1. It states in the other direction too that MiniSom's mexican hat is isotropic and so never had the separability bug this package shipped and fixed in 0.2.0. Every difference is presented as a convention choice, because that is what each one is. Two warnings a reader needs and would not otherwise get: `pip install sompy` installs an unrelated 2016 package by a different author, and `asv machine` writes the host's name, CPU model, RAM and kernel version to disk. No hostname appears in the recorded machine description. CPU model, core count, OS and Python version are enough to read the numbers against. The README gains a Benchmarks subsection pointing here. The CHANGELOG entry goes under Unreleased and opens by saying nothing in the shipped package changed, since a reader scanning for behaviour changes should be able to stop there. Also folds the one 105-character line left in the 0.6.1 entry, rather than spending a separate change on it. --- CHANGELOG.md | 42 +++- README.md | 15 ++ .../comparison-with-som-libraries.md | 188 ++++++++++++++++++ mkdocs.yml | 1 + 4 files changed, 244 insertions(+), 2 deletions(-) create mode 100644 docs/explanation/comparison-with-som-libraries.md diff --git a/CHANGELOG.md b/CHANGELOG.md index 6986d6f..a2611dd 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -5,6 +5,44 @@ All notable changes to this project are documented in this file. The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.1.0/), and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). +## [Unreleased] + +Nothing in the shipped package changed. This is benchmarking and the tests behind it. + +### Added + +- **A published comparison against MiniSom and SOMPY**, at + [Comparison with MiniSom and SOMPY](https://andremsouza.github.io/python-som/explanation/comparison-with-som-libraries/). + Both peers implement Kohonen's Eq. (8) and cite the same sources this package works from, so the + comparison is between three implementations of one published algorithm rather than between three + different algorithms. Getting there needed eight controls, because the packages look far more + interchangeable than they are: MiniSom's `train_batch` is stepwise training rather than batch, + only the gaussian is the same function in all three, and SOMPY's + `calculate_quantization_error` returns an elementwise mean absolute error rather than the usual + mean Euclidean distance. + + Measured on batch training, this package is 1.3x to 2.0x faster than SOMPY and 1.3x to 1.6x + **slower** than MiniSom, the latter growing with map size. `batch_update` walks every node in + Python, 108,000 `einsum` calls on a 60x60 map over 30 iterations, where MiniSom's node-side update + is a single vectorised divide. Stepwise training is within 10% of MiniSom's. + +- **`tests/test_minisom_agreement.py`**, which checks the neighborhood functions, the decay, the + best-matching-unit search and both training loops against MiniSom on every run. Trained models + agree to 2.8e-16 relative for stepwise and 1.3e-15 for batch, which is what makes timing the two + against each other meaningful. It also pins the two neighborhoods that deliberately *disagree*, so + neither is later "fixed" into the other. Same reasoning as + `tests/test_linalg_matches_sklearn.py`, which once found a real 5.8% defect in this package. + +- **An asv suite in `asv_benchmarks/`** tracking this package against its own git history, which is + what numpy, scipy, pandas and scikit-learn all use asv for. Sixteen benchmarks over the training + loops, the neighborhood kernel, the matching functions, the SVD and the artifact round trip. No + timing gates anything: CI executes each benchmark once and discards the numbers, because a shared + runner cannot measure anything. + +- **A `bench` extra** holding `asv` and `minisom`, kept out of `dev` because no gating job needs + them. SOMPY is deliberately absent and cannot be added: it uses `np.Inf`, removed in NumPy 2.0, at + class-definition time, so it gets an interpreter of its own that is built by hand. + ## [0.6.1] - 2026-07-30 ### Fixed @@ -20,8 +58,8 @@ adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). The upgrade notes were also short of two things worth knowing. 0.4.0 changes linear-initialization results for data far from the origin, where the previous PCA path was wrong by up to 5.8%, and it removed pandas and scikit-learn as runtime dependencies, which matters to - anyone who was importing either transitively. Both are now in the README, alongside a note that 0.5.0's - string deprecation was withdrawn so anyone who saw the warning knows to stop migrating. + anyone who was importing either transitively. Both are now in the README, alongside a note that + 0.5.0's string deprecation was withdrawn so anyone who saw the warning knows to stop migrating. ## [0.6.0] - 2026-07-30 diff --git a/README.md b/README.md index 17ce664..56b9dc4 100644 --- a/README.md +++ b/README.md @@ -134,6 +134,21 @@ pre-commit install # optional, run the gates on commit If you use the SonarQube for IDE (SonarLint) VS Code extension, it will also apply Sonar's Python rules locally; the ruff configuration is set up to cover most of the same ground. +### Benchmarks + +Hand-run, never part of the test suite: a timing assertion on shared hardware measures noise. + +```bash +uv run python benchmarks/bench_vs_minisom.py # vs MiniSom, agreement verified before timing +uv run python benchmarks/bench_batch.py # the neighborhood kernel against evaluating per node +cd asv_benchmarks && uv run --extra bench asv continuous master HEAD # this package across commits +``` + +Results and the method behind them are in +[Comparison with MiniSom and SOMPY](https://andremsouza.github.io/python-som/explanation/comparison-with-som-libraries/). +The SOMPY comparison needs an interpreter of its own, since SOMPY cannot be imported on NumPy 2; +`benchmarks/bench_vs_sompy.py` prints the setup command and installs nothing. + ## References Based on: diff --git a/docs/explanation/comparison-with-som-libraries.md b/docs/explanation/comparison-with-som-libraries.md new file mode 100644 index 0000000..aea08dd --- /dev/null +++ b/docs/explanation/comparison-with-som-libraries.md @@ -0,0 +1,188 @@ +# How this compares to MiniSom and SOMPY + +Three Python packages implement Kohonen's self-organizing map. This page measures python-som against +the other two, and spends most of its length on why that measurement is harder than it looks. + +The short version: on batch training python-som is **1.3x to 2.0x faster than SOMPY** and **1.3x to +1.6x slower than MiniSom**, with the gap against MiniSom growing as the map grows. On stepwise +training python-som and MiniSom are within about 10% of each other. SOMPY has no stepwise mode. + +## Why a naive comparison is wrong + +The three packages look interchangeable. python-som's method names deliberately mirror MiniSom's, +so `activate`, `winner`, `quantization`, `quantization_error` and `get_weights` all line up. +Underneath, some of them do different things. + +**MiniSom's `train_batch` is not batch training.** It is stepwise training in sequential sample +order, the counterpart of `mode="sequential"` here. Kohonen's Eq. (8) lives in +`train_batch_offline`. +Reaching for the name that matches would compare a weighted mean against per-sample gradient steps +and produce a number that looks plausible. + +**Only the gaussian is the same function in all three.** MiniSom writes it as an outer product of +two one-dimensional gaussians, which is valid because the exponential is separable; SOMPY applies +it to an already-squared grid distance. They agree to one unit in the last place. The other two +neighborhoods do not agree at all: + +| | python-som | MiniSom 2.3.6 | SOMPY @6aca604 | +| --- | --- | --- | --- | +| gaussian | $e^{-r^2/2\sigma^2}$ | separable product, same function | same function | +| mexican hat | $(1-u)e^{-u}$, zero at $\sqrt{2}\sigma$ | $(1-2u)e^{-u}$, zero at $\sigma$ | not implemented | +| bubble | Chebyshev, $\max(\|dx\|,\|dy\|) \le \mathrm{round}(\sigma)$ | Chebyshev, strict, unrounded | Euclidean disc | + +At $\sigma = 1$ the bubble selects 9 nodes here and 1 node in MiniSom. That is not a rounding +difference, and a benchmark using the bubble at equal $\sigma$ would compare a neighbourhood against +a point update. So every timing below uses the gaussian. + +None of these are defects. They are convention choices, and each package is internally consistent. +Worth stating plainly in one direction: **MiniSom's mexican hat is isotropic**, so it never had the +separability bug that python-som shipped in its first contributed version and fixed in 0.2.0. See +[Why isotropy matters](why-isotropy-matters.md) for what that bug was. + +## The controls + +With the differences known, the libraries can be driven into computing the same thing. Every control +below is asserted before any timing is reported, and the MiniSom ones are also checked on every CI +run by `tests/test_minisom_agreement.py`. + +1. **The same initial models are injected into all arms.** Seeding cannot do this: the three use + different generators, and their PCA initializations are three different functions. +2. **The gaussian only.** +3. **The same radius schedule.** `asymptotic_decay` is character-for-character identical in + python-som and MiniSom. SOMPY can only run `linspace(start, end, n)`, so python-som matches that + instead. +4. **Cases where python-som's radius floor never engages**, since neither peer has one. +5. **MiniSom's batch learning rate pinned to a constant 1.0.** `train_batch_offline` relaxes toward + the Eq. (8) mean by a decaying $\eta$ where Eq. (8) assigns it outright. +6. **Sequential rather than random sample order**, because `arange(T) % n` and + `np.resize(arange(n), T)` are the same sequence. The random modes differ genuinely and no seed + reconciles them. +7. **SOMPY's normalization turned off.** At its `'var'` default it would train on z-scored data. +8. **The quantization error recomputed under one definition for every arm.** This is not optional + for SOMPY: its `calculate_quantization_error` returns `mean(|x - m|)` over every feature, an + elementwise mean absolute error, where the standard definition is the mean Euclidean distance per + sample. SOMPY disagrees with itself here, since the value it logs during training *is* an L2 + distance. + +With those in place the trained models agree, which is what makes the timings mean anything: + +| against | agreement | why not exact | +| --- | --- | --- | +| MiniSom, stepwise | 2.8e-16 relative | floating-point association order | +| MiniSom, batch | 1.3e-15 relative | the two accumulate Eq. (8) in a different order | +| SOMPY, batch | 1.4e-07 relative | SOMPY rounds its codebook to six decimals every iteration | + +## The numbers + +Medians of 9 interleaved repeats. python-som 0.6.1, MiniSom 2.3.6 on NumPy 2.5.1, SOMPY @6aca604 on +NumPy 1.26.4, CPython 3.12.13, Linux, Intel Core Ultra 9 275HX. + +### Against MiniSom + +| map | samples | features | mode | python-som | MiniSom | result | +| --- | --- | --- | --- | --- | --- | --- | +| 20x20 | 200 | 4 | batch | 232.5 ms | 239.2 ms | 1.03x faster | +| 40x40 | 300 | 6 | batch | 1016.6 ms | 758.7 ms | **1.34x slower** | +| 60x60 | 400 | 8 | batch | 2885.2 ms | 1763.0 ms | **1.64x slower** | +| 20x20 | 200 | 4 | sequential | 1.2 ms | 1.4 ms | 1.12x faster | +| 40x40 | 300 | 6 | sequential | 2.4 ms | 2.5 ms | 1.06x faster | +| 60x60 | 400 | 8 | sequential | 4.4 ms | 4.7 ms | 1.08x faster | + +**python-som's batch training is slower, and the reason is structural rather than incidental.** +`batch_update` walks every node in Python: on a 60x60 map over 30 iterations that is 108,000 +`einsum` calls. MiniSom's node-side update is a single vectorised divide, and its Python loop runs +over the 400 samples instead. So the ratio tracks nodes against samples, which is why the gap is +absent at 20x20 and 1.64x at 60x60. Vectorising the node loop is the obvious fix, and has not been +done. + +### Against SOMPY + +| map | samples | features | python-som | SOMPY | result | +| --- | --- | --- | --- | --- | --- | +| 20x20 | 200 | 4 | 145.4 ms | 188.3 ms | 1.30x faster | +| 40x40 | 300 | 6 | 727.7 ms | 1264.6 ms | 1.74x faster | +| 60x60 | 400 | 8 | 1906.4 ms | 3796.6 ms | 1.99x faster | + +Batch only, since SOMPY implements nothing else: `SOMFactory.build` accepts `training='seq'` and +ignores it. + +Cases stop at 60x60 because SOMPY materialises the full node-by-node neighborhood, a `(3600, 3600)` +matrix at 104 MB. At 100x100 it would be 800 MB. python-som stores one `(2X-1, 2Y-1)` kernel +instead, 198 KB at 80x80, so its neighborhood memory grows linearly in node count where SOMPY's +grows quadratically. + +### Initialization + +The one default the three do not share, isolated by holding the algorithm and every hyperparameter +fixed and letting each library seed its own models. Quantization error after training, on data +scaled by 10 and offset 100 from the origin: + +| map | python-som | MiniSom | better | +| --- | --- | --- | --- | +| 20x20 | 1.6797 | 1.1157 | MiniSom | +| 40x40 | 0.2855 | 0.6265 | python-som | +| 60x60 | 0.0937 | 0.5009 | python-som | + +The three initializers place models on the plane of the first two principal components and differ in +how far apart. python-som uses $\bar{x} + c_1\sqrt{\lambda_1}v_1 + c_2\sqrt{\lambda_2}v_2$, scaling +each direction by the data's actual extent along it, from an SVD of the centred matrix. MiniSom uses +unit eigenvectors, so the initial sheet is one unit across whatever the data's spread. SOMPY scales +by $\lambda$ rather than $\sqrt{\lambda}$ and uses a randomized PCA, which makes its initialization +non-deterministic. + +## Reproducing this + +```bash +uv sync --extra dev --extra bench +uv run python benchmarks/bench_vs_minisom.py +``` + +SOMPY needs an interpreter of its own and you have to build it by hand, because it cannot be +installed alongside python-som at all: it uses `np.Inf` at class-definition time, and that was +removed in NumPy 2.0. `benchmarks/bench_vs_sompy.py` prints this command if the environment is +missing, and installs nothing itself. + +```bash +uv venv .venv-sompy --python 3.12 +uv pip install --python .venv-sompy \ + "numpy==1.26.4" "scipy==1.13.1" "scikit-learn==1.5.2" "scikit-image==0.24.0" \ + "numexpr==2.10.1" "joblib==1.4.2" "matplotlib==3.9.2" \ + "SOMPY @ git+https://github.com/sevamoo/SOMPY@6aca604b06e5eea1391ecf507810c7aabafc3f8b" +uv run python benchmarks/bench_vs_sompy.py +``` + +!!! warning "`pip install sompy` gets a different package" + + The SOMPY compared here is [sevamoo/SOMPY](https://github.com/sevamoo/SOMPY), which was never + published to PyPI. The PyPI project named `sompy` is version 0.1.1 by a different author, two + releases in 2016, and is unrelated. The pin above is a commit SHA because the repository has no + tags. + +To track python-som against its own history rather than against a peer, which is a different +question and wants a different tool: + +```bash +cd asv_benchmarks +uv run --extra bench asv machine --yes +uv run --extra bench asv continuous master HEAD +``` + +!!! note "`asv machine` records your hardware" + + It writes your host's name, CPU model, core count, RAM and kernel version to + `.asv/results//machine.json`. That directory is gitignored. + +## Caveats + +Timings come from one machine and are not portable; the ratios are more durable than the +milliseconds, and both are worth re-measuring rather than quoting. Nothing here is a CI gate, +because shared runners cannot measure anything: an earlier version of this comparison swung from +2.56x to 1.01x between two runs on an otherwise idle laptop. + +The peers are pinned, and MiniSom's development branch has already moved past its release with a +numba-compiled batch path that would change these numbers substantially. + +Both peers are worth using. MiniSom is faster at batch training on larger maps and has hexagonal +topologies, which python-som does not. SOMPY has clustering and visualization helpers built in. What +python-som offers against them is NumPy as its only runtime dependency, a scikit-learn estimator +interface, pickle-free artifacts with provenance, and type information. diff --git a/mkdocs.yml b/mkdocs.yml index 7767466..a043061 100644 --- a/mkdocs.yml +++ b/mkdocs.yml @@ -45,6 +45,7 @@ nav: - Why isotropy matters: explanation/why-isotropy-matters.md - Batch vs stepwise: explanation/batch-vs-stepwise.md - Artifact safety: explanation/artifact-safety.md + - Comparison with MiniSom and SOMPY: explanation/comparison-with-som-libraries.md plugins: - search