1. مقدمه
این پیوست، قراردادهای علمی فصل اصلی را به یک جریان اجرایی قابل بازتولید تبدیل میکند. هدف، تکرار مبانی نظری نیست؛ تمام تصمیمهای کدی به نمادگذاری و محدودیتهای تثبیتشده در فصل ارجاع دارند. در پیادهسازی حاضر، بردار ورودی با x، بردار کد واحد j باm، مختصات شبکهای با rj، واحد با بهترین تطابق با c(x)، تابع همسایگی با hcj، نرخ یادگیری با α(t) و شعاع همسایگی با σ(t) نمایش داده میشوند.
دو قرارداد ارزیابی فصل بدون تغییر حفظ شدهاند:
در این رابطه، QE میانگین فاصله اقلیدسی بدون توان دوم تا BMU است؛ بنابراین نباید بدون اعلام تغییر قرارداد با MSE یا مجموع فاصلههای مربعشده جایگزین شود.
حالا در این رابطه، Adj(c1n,c2n ) نشان میدهد نخستین و دومین BMU در شبکه مجاورند. در کد پایه، مجاورت شبکه مستطیلی بهصورت چهارهمسایگیِ مشترک در ضلع تعریف شده است. پس در توپولوژی حلقوی، لبههای مقابل نیز همسایهاند. این قرارداد همراه هر نتیجه ثبت میشود تا TE میان کتابخانهها بهاشتباه مقایسه نشود.
بسته همراه پیوست شامل کد منبع، آزمونها، خروجیهای CSV و JSON، تصاویر، نسخه دقیق محیط، seedها و SHA-256 همه خروجیهاست. برچسبهای کلاس در هیچیک از مطالعات برای آموزش SOM استفاده نشدهاند؛ برچسبها فقط در ارزیابی پسینی به کار میروند.
.
2. پیشنیازهای عملی و ساختار داده ورودی
2.1 محیط اجرایی
دو دسته نسخه ثبت شدهاند. «نسخههای پایدار بررسیشده» برای نصب جدید پیشنهاد میشوند، ولی اعداد گزارششده این پیوست با «محیط دقیق آزمونشده» تولید شدهاند.
| مؤلفه | نسخه پایدار بررسیشده در ۵ اوت ۲۰۲۶ | نسخه دقیق تولید نتایج |
| Python | 3.14.6 | 3.13.5 |
| NumPy | 2.5.1 | 2.3.5 |
| SciPy | 1.18.0 | 1.17.0 |
| scikit-learn | 1.9.0 | 1.8.0 |
| pandas | 3.0.4 | 2.2.3 |
| Matplotlib | 3.11.0 | 3.10.8 |
| pytest | 9.1.1 | 9.0.2 |
نسخههای پایدار از صفحات انتشار رسمی پروژهها استخراج شدهاند. برای بازتولید عددی نتایج کتاب باید از requirements-tested.txt استفاده شود؛ ارتقا به نسخههای جدید ممکن است بهدلیل تفاوتهای محاسبات ممیز شناور، مقداردهی اولیه PCA یا پیادهسازی K-means در چند رقم اعشار تفاوت ایجاد کند.
نصب محیط آزمونشده
python -m venv .venv
source .venv/bin/activate # Windows: .venv\Scripts\activate
python -m pip install -r requirements-tested.txt
python -m pip install -e .
python run_all.py
اجرای run_all.py سه مطالعه را اجرا میکند، تصاویر و جدولها را میسازد، آزمونها را اجرا میکند، نسخهها و seedها را در environment.json ثبت میکند و برای خروجیها checksum میسازد.
2.2 شمای ورودی
ورودی پیادهسازی پایه آرایهای عددی و دوبعدی با شکل زیر است:
- هر سطر یک مشاهده و هر ستون یک ویژگی است.
- داده باید float64 یا قابل تبدیل به آن باشد.
- مقدار NaN یا بینهایت مستقیماً پذیرفته نمیشود.
- پیشپردازش فقط روی مجموعه آموزش برازش و سپس بدون بازبرازش روی آزمون اعمال میشود.
- برچسب، شناسه رکورد، زمان و ستونهای غیرعددی نباید بهصورت خام وارد فاصله اقلیدسی شوند.
.
2.3 پیشپردازش ضروری
مقیاس ویژگیها
SOM پایه از فاصله اقلیدسی استفاده میکند. اگر یک ویژگی دامنهای چند هزار برابر ویژگی دیگر داشته باشد، BMU تقریباً بر اساس همان ویژگی تعیین میشود. گزینه پیشفرض مناسب برای دادههای معمول StandardScaler است. برای توزیعهای سنگیندُم یا داده پرت، RobustScaler و برش صدکیِ توجیهشده میتوانند مناسبتر باشند. این مداخله باید در گزارش نتیجه ثبت شود، زیرا فضای فاصله را تغییر میدهد.
داده ناقص
جانشانی باید پیش از SOM انجام شود. در مطالعات این پیوست از میانه آموزش استفاده شده است. محاسبه فاصله روی تعداد متغیر ویژگیهای مشاهدهشده در پیادهسازی پایه پشتیبانی نمیشود؛ چنین روشی نیازمند تعریف نرمالسازی و ارزیابی مستقل است.
ویژگیهای همبسته و ابعاد بالا
PCA میتواند هزینه و اثر ویژگیهای همبسته را کاهش دهد، ولی مؤلفههای بیشینهواریانس الزاماً اطلاعات کلاس اقلیت یا ساختار نادر را حفظ نمیکنند. مطالعه دوم دقیقاً یک حالت شکست را نشان میدهد که در آن PCA مقدار را کاهش میدهد، اما بازخوانی کلاس اقلیت را بهشدت تضعیف میکند.
.
2.4 قراردادهای بازتولیدپذیری
- تمام seedها ثابت و در environment.json ثبت شدهاند.
- مقداردهی اولیه PCA با قرارداد علامت قطعی SVD پیادهسازی شده است تا جابهجایی ترتیب سطرها نقشه را صرفاً بهدلیل ابهام علامت مؤلفه منعکس نکند.
- در تساوی فاصله BMU، واحد دارای اندیس کوچکتر انتخاب میشود.
- آموزش دستهای داده را در بلوکها پردازش میکند و ماتریس کامل را در حافظه نگه نمیدارد.
- برای نسخه آنلاین، ترتیب نمونهها با مولد NumPy و seed مشخص کنترل میشود.
- مقایسه چند پیشپردازش بر پایه تنها وقتی معتبر است که فضای ویژگی و مقیاس یکسان باشد.
.
3. پیادهسازی پایه از صفر
3.1 تصمیمهای طراحی
کلاس SelfOrganizingMap فقط به NumPy وابسته است و دو شیوه آموزش را پشتیبانی میکند:
- training=”online”: اجرای قانون کلاسیک نمونهبهنمونه با و ؛
- training=”batch”: تعیین BMU همه نمونهها و محاسبه میانگین وزندار، بدون نرخ یادگیری صریح.
دو توپولوژی در نسخه پایه وجود دارد: شبکه مستطیلی صفحهای و شبکه حلقوی. تابع همسایگی میتواند گاوسی یا حبابی باشد. مقداردهی اولیه PCA، نمونهمحور یا یکنواخت پشتیبانی شده است. فاصله دادهای در این نسخه عمداً اقلیدسی است؛ جایگزینی DTW یا عدمشباهت گرافی بدون تعریف prototype سازگار، یک امتداد جداگانه محسوب میشود.
.
3.2 کد کامل ماژول پایه
کد زیر همان فایل src/som_operational/som.py در بسته بازتولیدپذیری است. کنترل ورودی، آموزش آنلاین و دستهای، محاسبه و ، U-matrix، hit map، صفحههای مؤلفه و برچسبگذاری پسینی واحدها را شامل میشود.
"""A reproducible NumPy implementation of the Self-Organizing Map (SOM).
The implementation follows the notation used in the accompanying chapter:
* input vector: x
* prototype/codebook vector of unit j: m_j
* lattice coordinate of unit j: r_j
* best matching unit: c(x)
* neighborhood function: h_cj
* online learning rate: alpha(t)
* neighborhood radius: sigma(t)
* quantization error: QE, mean *non-squared* Euclidean distance
* topographic error: TE, based on first/second BMU adjacency
Only complete, finite, numeric input is accepted. Missing-value handling and
feature scaling belong to the preprocessing pipeline and are intentionally not
hidden inside the estimator.
"""
from __future__ import annotations
from dataclasses import dataclass
from typing import Literal, Sequence
import numpy as np
from numpy.typing import ArrayLike, NDArray
TrainingMode = Literal["online", "batch"]
Topology = Literal["rectangular", "toroidal"]
Neighborhood = Literal["gaussian", "bubble"]
Initialization = Literal["pca", "random_samples", "random_uniform"]
@dataclass(frozen=True)
class SOMHistoryEntry:
"""Metrics recorded after an epoch."""
epoch: int
quantization_error: float
topographic_error: float
max_codebook_shift: float
learning_rate: float | None
sigma: float
class SelfOrganizingMap:
"""Self-Organizing Map with online and batch learning.
Parameters
----------
map_shape:
``(rows, columns)`` of the two-dimensional lattice.
n_epochs:
Maximum number of complete training epochs.
training:
``"online"`` for sample-wise Kohonen updates or ``"batch"`` for
weighted batch updates.
learning_rate:
``(alpha_start, alpha_end)`` for online training. It is recorded but
not used by batch training, which has no explicit learning-rate term.
sigma:
``(sigma_start, sigma_end)``. If the first value is ``None``, half the
largest map dimension is used.
topology:
``"rectangular"`` for an ordinary planar map or ``"toroidal"`` to
connect opposite boundaries.
neighborhood:
Gaussian or bubble neighborhood on the lattice.
initialization:
PCA-linear, sampled observations, or uniform within feature bounds.
shuffle:
Shuffle row order at each online epoch. It has no mathematical effect
on an exact batch epoch and is therefore ignored in batch mode.
random_state:
Seed for initialization and online sample order.
tol, patience:
Optional early stopping on maximum prototype movement.
batch_size:
Chunk size used for BMU search and batch sufficient-statistic
accumulation. This avoids materializing an ``N x M`` distance matrix.
record_history:
Compute QE and TE after every epoch. Disable for very large data.
Notes
-----
* Rectangular-grid adjacency used by :meth:`topographic_error` is the
four-neighbour (edge-sharing) convention. Toroidal adjacency uses the
same convention with wrapped boundaries.
* The codebooks are ordinary vectors in Euclidean feature space. Replacing
Euclidean distance with a structured dissimilarity such as DTW requires
a compatible prototype representation/update and is outside this base
implementation.
"""
def __init__(
self,
map_shape: tuple[int, int] = (8, 6),
*,
n_epochs: int = 100,
training: TrainingMode = "batch",
learning_rate: tuple[float, float] = (0.5, 0.05),
sigma: tuple[float | None, float] = (None, 0.8),
topology: Topology = "rectangular",
neighborhood: Neighborhood = "gaussian",
initialization: Initialization = "pca",
shuffle: bool = True,
random_state: int | None = 0,
tol: float = 1e-5,
patience: int = 5,
batch_size: int = 512,
record_history: bool = True,
) -> None:
self.map_shape = map_shape
self.n_epochs = n_epochs
self.training = training
self.learning_rate = learning_rate
self.sigma = sigma
self.topology = topology
self.neighborhood = neighborhood
self.initialization = initialization
self.shuffle = shuffle
self.random_state = random_state
self.tol = tol
self.patience = patience
self.batch_size = batch_size
self.record_history = record_history
self.codebook_: NDArray[np.float64] | None = None
self.lattice_: NDArray[np.float64] | None = None
self.history_: list[SOMHistoryEntry] = []
self.n_features_in_: int | None = None
self.n_iter_: int = 0
self._validate_hyperparameters()
@property
def n_units(self) -> int:
return int(self.map_shape[0] * self.map_shape[1])
def _validate_hyperparameters(self) -> None:
rows, cols = self.map_shape
if not isinstance(rows, int) or not isinstance(cols, int) or rows < 1 or cols < 1:
raise ValueError("map_shape must contain two positive integers.")
if not isinstance(self.n_epochs, int) or self.n_epochs < 1:
raise ValueError("n_epochs must be a positive integer.")
if self.training not in {"online", "batch"}:
raise ValueError("training must be 'online' or 'batch'.")
if self.topology not in {"rectangular", "toroidal"}:
raise ValueError("topology must be 'rectangular' or 'toroidal'.")
if self.neighborhood not in {"gaussian", "bubble"}:
raise ValueError("neighborhood must be 'gaussian' or 'bubble'.")
if self.initialization not in {"pca", "random_samples", "random_uniform"}:
raise ValueError("Unsupported initialization method.")
alpha0, alphaf = self.learning_rate
if alpha0 <= 0 or alphaf <= 0 or alphaf > alpha0:
raise ValueError("learning_rate must be positive and non-increasing.")
sigma0, sigmaf = self.sigma
if sigma0 is not None and sigma0 <= 0:
raise ValueError("Initial sigma must be positive or None.")
if sigmaf <= 0:
raise ValueError("Final sigma must be positive.")
if sigma0 is not None and sigmaf > sigma0:
raise ValueError("sigma must be non-increasing.")
if self.tol < 0:
raise ValueError("tol must be non-negative.")
if not isinstance(self.patience, int) or self.patience < 1:
raise ValueError("patience must be a positive integer.")
if not isinstance(self.batch_size, int) or self.batch_size < 1:
raise ValueError("batch_size must be a positive integer.")
@staticmethod
def _as_2d_finite(X: ArrayLike) -> NDArray[np.float64]:
arr = np.asarray(X, dtype=np.float64)
if arr.ndim != 2:
raise ValueError("X must be a two-dimensional numeric array.")
if arr.shape[0] < 1 or arr.shape[1] < 1:
raise ValueError("X must contain at least one row and one feature.")
if not np.all(np.isfinite(arr)):
raise ValueError("X contains NaN or infinite values; preprocess them explicitly.")
return np.ascontiguousarray(arr)
def _require_fitted(self) -> tuple[NDArray[np.float64], NDArray[np.float64]]:
if self.codebook_ is None or self.lattice_ is None:
raise RuntimeError("The SOM is not fitted.")
return self.codebook_, self.lattice_
def _make_lattice(self) -> NDArray[np.float64]:
rows, cols = self.map_shape
return np.array([(r, c) for r in range(rows) for c in range(cols)], dtype=np.float64)
def _initialize_codebook(
self, X: NDArray[np.float64], rng: np.random.Generator
) -> NDArray[np.float64]:
M, d = self.n_units, X.shape[1]
if self.initialization == "random_samples":
idx = rng.choice(X.shape[0], size=M, replace=X.shape[0] < M)
return X[idx].copy()
if self.initialization == "random_uniform":
lo, hi = X.min(axis=0), X.max(axis=0)
return rng.uniform(lo, hi, size=(M, d))
# PCA-linear initialization. SVD is stable and works for d=1.
mean = X.mean(axis=0)
centered = X - mean
_, s, vt = np.linalg.svd(centered, full_matrices=False)
# Resolve the arbitrary SVD sign so row permutations do not mirror the map.
for k in range(vt.shape[0]):
pivot = int(np.argmax(np.abs(vt[k])))
if vt[k, pivot] < 0:
vt[k] *= -1.0
scale = s / np.sqrt(max(1, X.shape[0] - 1))
rows, cols = self.map_shape
row_axis = np.linspace(-1.0, 1.0, rows)
col_axis = np.linspace(-1.0, 1.0, cols)
codebook = np.tile(mean, (M, 1))
for unit, (r, c) in enumerate(self._make_lattice().astype(int)):
if vt.shape[0] >= 1:
coord1 = col_axis[c] if cols > 1 else row_axis[r]
codebook[unit] += coord1 * scale[0] * vt[0]
if vt.shape[0] >= 2:
coord2 = row_axis[r]
codebook[unit] += coord2 * scale[1] * vt[1]
return np.clip(codebook, X.min(axis=0), X.max(axis=0))
@staticmethod
def _geometric_schedule(start: float, end: float, progress: float) -> float:
if np.isclose(start, end):
return float(start)
return float(start * (end / start) ** progress)
def _epoch_parameters(self, epoch: int) -> tuple[float, float]:
progress = 1.0 if self.n_epochs == 1 else epoch / (self.n_epochs - 1)
alpha = self._geometric_schedule(
float(self.learning_rate[0]), float(self.learning_rate[1]), progress
)
sigma_start = (
max(self.map_shape) / 2.0 if self.sigma[0] is None else float(self.sigma[0])
)
sigma = self._geometric_schedule(sigma_start, float(self.sigma[1]), progress)
return alpha, sigma
def _lattice_squared_distances(self, bmu_indices: NDArray[np.int64]) -> NDArray[np.float64]:
_, lattice = self._require_fitted()
origins = lattice[bmu_indices]
delta = np.abs(origins[:, None, :] - lattice[None, :, :])
if self.topology == "toroidal":
shape = np.asarray(self.map_shape, dtype=np.float64)
delta = np.minimum(delta, shape[None, None, :] - delta)
return np.sum(delta * delta, axis=2)
def _neighborhood_weights(
self, bmu_indices: NDArray[np.int64], sigma: float
) -> NDArray[np.float64]:
squared = self._lattice_squared_distances(bmu_indices)
if self.neighborhood == "gaussian":
return np.exp(-squared / (2.0 * sigma * sigma))
return (squared <= sigma * sigma).astype(np.float64)
@staticmethod
def _squared_euclidean(
X: NDArray[np.float64], codebook: NDArray[np.float64]
) -> NDArray[np.float64]:
# max(..., 0) removes tiny negative values caused by roundoff.
return np.maximum(
np.sum(X * X, axis=1, keepdims=True)
+ np.sum(codebook * codebook, axis=1)[None, :]
- 2.0 * X @ codebook.T,
0.0,
)
def _bmu_indices_in_chunks(
self, X: NDArray[np.float64], *, second: bool = False
) -> NDArray[np.int64] | tuple[NDArray[np.int64], NDArray[np.int64]]:
codebook, _ = self._require_fitted()
first_parts: list[NDArray[np.int64]] = []
second_parts: list[NDArray[np.int64]] = []
for start in range(0, X.shape[0], self.batch_size):
chunk = X[start : start + self.batch_size]
dist2 = self._squared_euclidean(chunk, codebook)
if second:
# Stable sort gives deterministic tie handling by unit index.
order = np.argsort(dist2, axis=1, kind="stable")[:, :2]
first_parts.append(order[:, 0].astype(np.int64))
second_parts.append(order[:, 1].astype(np.int64))
else:
first_parts.append(np.argmin(dist2, axis=1).astype(np.int64))
first_arr = np.concatenate(first_parts)
if not second:
return first_arr
return first_arr, np.concatenate(second_parts)
def fit(self, X: ArrayLike) -> "SelfOrganizingMap":
"""Fit codebook vectors to complete numeric data."""
data = self._as_2d_finite(X)
self.n_features_in_ = data.shape[1]
self.lattice_ = self._make_lattice()
rng = np.random.default_rng(self.random_state)
self.codebook_ = self._initialize_codebook(data, rng)
self.history_ = []
stale_epochs = 0
for epoch in range(self.n_epochs):
alpha, sigma = self._epoch_parameters(epoch)
before = self.codebook_.copy()
if self.training == "online":
order = np.arange(data.shape[0])
if self.shuffle:
rng.shuffle(order)
for idx in order:
x = data[idx]
dist2 = np.sum((self.codebook_ - x) ** 2, axis=1)
c = int(np.argmin(dist2))
h = self._neighborhood_weights(np.array([c], dtype=np.int64), sigma)[0]
self.codebook_ += alpha * h[:, None] * (x - self.codebook_)
alpha_record: float | None = alpha
else:
numerator = np.zeros_like(self.codebook_)
denominator = np.zeros(self.n_units, dtype=np.float64)
for start in range(0, data.shape[0], self.batch_size):
chunk = data[start : start + self.batch_size]
dist2 = self._squared_euclidean(chunk, self.codebook_)
bmu = np.argmin(dist2, axis=1).astype(np.int64)
h = self._neighborhood_weights(bmu, sigma)
numerator += h.T @ chunk
denominator += h.sum(axis=0)
occupied = denominator > np.finfo(np.float64).eps
self.codebook_[occupied] = numerator[occupied] / denominator[occupied, None]
alpha_record = None
shift = float(np.max(np.linalg.norm(self.codebook_ - before, axis=1)))
if self.record_history:
self.history_.append(
SOMHistoryEntry(
epoch=epoch + 1,
quantization_error=self.quantization_error(data),
topographic_error=self.topographic_error(data),
max_codebook_shift=shift,
learning_rate=alpha_record,
sigma=sigma,
)
)
self.n_iter_ = epoch + 1
if shift <= self.tol:
stale_epochs += 1
if stale_epochs >= self.patience:
break
else:
stale_epochs = 0
return self
def predict_bmu(self, X: ArrayLike) -> NDArray[np.int64]:
"""Return flat BMU indices for rows of ``X``."""
data = self._as_2d_finite(X)
if self.n_features_in_ is not None and data.shape[1] != self.n_features_in_:
raise ValueError("X has a different number of features from the fitted data.")
return self._bmu_indices_in_chunks(data) # type: ignore[return-value]
def predict_coordinates(self, X: ArrayLike) -> NDArray[np.int64]:
"""Return ``(row, column)`` coordinates of each BMU."""
_, lattice = self._require_fitted()
return lattice[self.predict_bmu(X)].astype(np.int64)
def transform(self, X: ArrayLike) -> NDArray[np.float64]:
"""Return Euclidean distance from each row to every codebook vector."""
data = self._as_2d_finite(X)
codebook, _ = self._require_fitted()
if self.n_features_in_ is not None and data.shape[1] != self.n_features_in_:
raise ValueError("X has a different number of features from the fitted data.")
chunks = []
for start in range(0, data.shape[0], self.batch_size):
dist2 = self._squared_euclidean(data[start : start + self.batch_size], codebook)
chunks.append(np.sqrt(dist2))
return np.vstack(chunks)
def quantization_distances(self, X: ArrayLike) -> NDArray[np.float64]:
"""Distance of each observation to its BMU."""
distances = self.transform(X)
return np.min(distances, axis=1)
def quantization_error(self, X: ArrayLike) -> float:
"""Mean non-squared Euclidean distance to the BMU, matching the chapter."""
return float(np.mean(self.quantization_distances(X)))
def _are_adjacent(self, first: NDArray[np.int64], second: NDArray[np.int64]) -> NDArray[np.bool_]:
_, lattice = self._require_fitted()
a, b = lattice[first], lattice[second]
delta = np.abs(a - b)
if self.topology == "toroidal":
shape = np.asarray(self.map_shape, dtype=np.float64)
delta = np.minimum(delta, shape[None, :] - delta)
# Four-neighbour rectangular convention: exactly one unit along one axis.
return np.sum(delta, axis=1) == 1.0
def topographic_error(self, X: ArrayLike) -> float:
"""Fraction whose first and second BMUs are not lattice neighbours."""
data = self._as_2d_finite(X)
if self.n_units < 2:
return 0.0
first, second = self._bmu_indices_in_chunks(data, second=True) # type: ignore[misc]
return float(np.mean(~self._are_adjacent(first, second)))
def hit_map(self, X: ArrayLike) -> NDArray[np.int64]:
"""Count observations assigned to each unit, shaped like the map."""
hits = np.bincount(self.predict_bmu(X), minlength=self.n_units)
return hits.reshape(self.map_shape)
def u_matrix(self) -> NDArray[np.float64]:
"""Mean codebook distance to edge-sharing lattice neighbours."""
codebook, lattice = self._require_fitted()
result = np.zeros(self.n_units, dtype=np.float64)
counts = np.zeros(self.n_units, dtype=np.int64)
rows, cols = self.map_shape
for idx, (r_float, c_float) in enumerate(lattice):
r, c = int(r_float), int(c_float)
neighbours: set[tuple[int, int]] = set()
for dr, dc in ((-1, 0), (1, 0), (0, -1), (0, 1)):
rr, cc = r + dr, c + dc
if self.topology == "toroidal":
rr %= rows
cc %= cols
if 0 <= rr < rows and 0 <= cc < cols and (rr, cc) != (r, c):
neighbours.add((rr, cc))
for rr, cc in neighbours:
j = rr * cols + cc
result[idx] += np.linalg.norm(codebook[idx] - codebook[j])
counts[idx] += 1
result /= np.maximum(counts, 1)
return result.reshape(self.map_shape)
def component_planes(self) -> NDArray[np.float64]:
"""Return an array of shape ``(n_features, rows, columns)``."""
codebook, _ = self._require_fitted()
return codebook.T.reshape((codebook.shape[1], *self.map_shape))
def get_history_array(self) -> NDArray[np.float64]:
"""History as numeric columns: epoch, QE, TE, shift, alpha, sigma."""
if not self.history_:
return np.empty((0, 6), dtype=np.float64)
return np.asarray(
[
[
h.epoch,
h.quantization_error,
h.topographic_error,
h.max_codebook_shift,
np.nan if h.learning_rate is None else h.learning_rate,
h.sigma,
]
for h in self.history_
],
dtype=np.float64,
)
def majority_unit_labels(
bmu_indices: ArrayLike,
y: ArrayLike,
n_units: int,
*,
codebook: ArrayLike | None = None,
) -> NDArray[np.int64]:
"""Assign post-hoc majority labels to SOM units.
Empty units are filled from the nearest occupied codebook vector when a
codebook is supplied. Labels are used only *after* unsupervised training.
"""
bmu = np.asarray(bmu_indices, dtype=np.int64)
labels = np.asarray(y, dtype=np.int64)
if bmu.ndim != 1 or labels.ndim != 1 or bmu.shape[0] != labels.shape[0]:
raise ValueError("bmu_indices and y must be one-dimensional and aligned.")
if n_units < 1:
raise ValueError("n_units must be positive.")
unit_labels = np.full(n_units, -1, dtype=np.int64)
occupied: list[int] = []
for unit in range(n_units):
vals = labels[bmu == unit]
if vals.size:
counts = np.bincount(vals)
unit_labels[unit] = int(np.flatnonzero(counts == counts.max())[0])
occupied.append(unit)
if not occupied:
raise ValueError("No occupied units were found.")
if np.any(unit_labels < 0):
if codebook is None:
global_counts = np.bincount(labels)
unit_labels[unit_labels < 0] = int(np.argmax(global_counts))
else:
cb = np.asarray(codebook, dtype=np.float64)
occupied_arr = np.asarray(occupied, dtype=np.int64)
for unit in np.flatnonzero(unit_labels < 0):
nearest = occupied_arr[
np.argmin(np.sum((cb[occupied_arr] - cb[unit]) ** 2, axis=1))
]
unit_labels[unit] = unit_labels[nearest]
return unit_labels
def unit_cluster_assignments(
codebook: ArrayLike,
n_clusters: int,
*,
random_state: int = 0,
sample_weight: ArrayLike | None = None,
) -> NDArray[np.int64]:
"""Cluster codebooks for an explicit downstream partition.
``sample_weight`` can be the training hit count of each unit, which prevents
empty units from dominating the secondary k-means partition.
"""
from sklearn.cluster import KMeans
cb = np.asarray(codebook, dtype=np.float64)
if cb.ndim != 2:
raise ValueError("codebook must be a two-dimensional array.")
weights = None if sample_weight is None else np.asarray(sample_weight, dtype=np.float64)
if weights is not None and weights.shape != (cb.shape[0],):
raise ValueError("sample_weight must contain one value per codebook vector.")
model = KMeans(n_clusters=n_clusters, n_init="auto", random_state=random_state)
model.fit(cb, sample_weight=weights)
return model.labels_.astype(np.int64)
3.3.آزمونهای واحد
آزمونها پنج قرارداد علمی حساس را مستقیم کنترل میکنند:
- محاسبه BMU؛
- بازتولید دقیق مثال عددی فصل؛
- استفاده از فاصله بدون توان دوم در ؛
- تعریف چهارهمسایگی و wrap در ؛
- ناوردایی یک دوره آموزش دستهای نسبت به جایگشت سطرها.
همچنین رد داده ناقص، شکل U-matrix، بازتولیدپذیری seed و برچسبگذاری واحدهای خالی آزمون شدهاند. نتیجه محیط ثبتشده:
.......... [100%]
10 passed in 0.11s
کد کامل آزمونها:
from __future__ import annotations
import numpy as np
import pytest
from som_operational import SelfOrganizingMap, majority_unit_labels
def test_known_bmu() -> None:
som = SelfOrganizingMap(map_shape=(1, 2), n_epochs=1)
som.lattice_ = som._make_lattice()
som.codebook_ = np.array([[0.0, 0.0], [1.0, 1.0]])
som.n_features_in_ = 2
assert som.predict_bmu([[0.8, 0.9]]).tolist() == [1]
def test_chapter_single_update() -> None:
som = SelfOrganizingMap(
map_shape=(2, 2),
n_epochs=1,
training="online",
learning_rate=(0.5, 0.5),
sigma=(1.0, 1.0),
initialization="random_uniform",
shuffle=False,
tol=0.0,
patience=1,
record_history=False,
)
som.lattice_ = som._make_lattice()
som.codebook_ = np.array([[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]])
som.n_features_in_ = 2
# Execute the same one-step update as the chapter without reinitialization.
x = np.array([0.8, 0.2])
c = int(np.argmin(np.sum((som.codebook_ - x) ** 2, axis=1)))
h = som._neighborhood_weights(np.array([c]), 1.0)[0]
som.codebook_ += 0.5 * h[:, None] * (x - som.codebook_)
expected = np.array(
[
[0.24261226, 0.06065307],
[0.9, 0.1],
[0.14715178, 0.85284822],
[0.93934693, 0.75738774],
]
)
np.testing.assert_allclose(som.codebook_, expected, rtol=1e-6, atol=1e-6)
def test_quantization_error_is_non_squared_distance() -> None:
som = SelfOrganizingMap(map_shape=(1, 1), n_epochs=1)
som.lattice_ = som._make_lattice()
som.codebook_ = np.array([[0.0, 0.0]])
som.n_features_in_ = 2
assert som.quantization_error([[3.0, 4.0], [0.0, 0.0]]) == pytest.approx(2.5)
def test_rectangular_topographic_error_uses_four_neighbours() -> None:
som = SelfOrganizingMap(map_shape=(2, 2), n_epochs=1, topology="rectangular")
som.lattice_ = som._make_lattice()
som.codebook_ = np.array([[0.0], [10.0], [9.0], [1.0]])
som.n_features_in_ = 1
# x=0.4 has first BMU unit 0 and second BMU unit 3; diagonal -> not adjacent.
assert som.topographic_error([[0.4]]) == pytest.approx(1.0)
def test_toroidal_boundaries_are_adjacent() -> None:
som = SelfOrganizingMap(map_shape=(1, 4), n_epochs=1, topology="toroidal")
som.lattice_ = som._make_lattice()
som.codebook_ = np.array([[0.0], [10.0], [9.0], [1.0]])
som.n_features_in_ = 1
# Units 0 and 3 are adjacent only under toroidal wrapping.
assert som.topographic_error([[0.4]]) == pytest.approx(0.0)
def test_batch_training_is_row_order_invariant() -> None:
X = np.array([[-1.0, -1.0], [-0.8, -1.2], [1.0, 1.0], [0.9, 1.1]])
common = dict(
map_shape=(2, 2),
n_epochs=8,
training="batch",
sigma=(1.5, 0.7),
initialization="pca",
random_state=7,
tol=0.0,
record_history=False,
)
first = SelfOrganizingMap(**common).fit(X)
second = SelfOrganizingMap(**common).fit(X[[2, 0, 3, 1]])
np.testing.assert_allclose(first.codebook_, second.codebook_, atol=1e-12)
def test_reproducibility_with_same_seed() -> None:
rng = np.random.default_rng(123)
X = rng.normal(size=(40, 3))
kwargs = dict(
map_shape=(3, 3),
n_epochs=5,
training="online",
initialization="random_samples",
random_state=99,
tol=0.0,
record_history=False,
)
a = SelfOrganizingMap(**kwargs).fit(X)
b = SelfOrganizingMap(**kwargs).fit(X)
np.testing.assert_allclose(a.codebook_, b.codebook_)
def test_nan_is_rejected() -> None:
with pytest.raises(ValueError, match="NaN"):
SelfOrganizingMap().fit([[1.0, np.nan]])
def test_u_matrix_and_component_plane_shapes() -> None:
X = np.array([[0.0, 0.0], [1.0, 1.0], [0.0, 1.0], [1.0, 0.0]])
som = SelfOrganizingMap(map_shape=(2, 3), n_epochs=3, random_state=2).fit(X)
assert som.u_matrix().shape == (2, 3)
assert som.component_planes().shape == (2, 2, 3)
assert som.hit_map(X).sum() == X.shape[0]
def test_majority_labels_fill_empty_units() -> None:
bmu = np.array([0, 0, 2, 2])
y = np.array([1, 1, 0, 0])
codebook = np.array([[0.0], [0.1], [1.0]])
labels = majority_unit_labels(bmu, y, 3, codebook=codebook)
assert labels.tolist() == [1, 1, 0]
4. مثال آموزشی کوچک
این مثال دقیقاً با مثال عددی بخش ۹ فصل اصلی همتراز است. شبکه ، ورودی ، نرخ یادگیری و شعاع دارد.
4.1 کد کامل
"""One exact online SOM update matching the numerical example in the chapter."""
from __future__ import annotations
from pathlib import Path
import sys
import numpy as np
import pandas as pd
ROOT = Path(__file__).resolve().parents[1]
sys.path.insert(0, str(ROOT / "src"))
from som_operational import SelfOrganizingMap # noqa: E402
def run(output_dir: Path) -> pd.DataFrame:
output_dir.mkdir(parents=True, exist_ok=True)
initial = np.array(
[[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]], dtype=float
)
x = np.array([0.8, 0.2], dtype=float)
alpha = 0.5
sigma = 1.0
som = SelfOrganizingMap(
map_shape=(2, 2),
n_epochs=1,
training="online",
learning_rate=(alpha, alpha),
sigma=(sigma, sigma),
initialization="random_uniform",
shuffle=False,
record_history=False,
)
som.lattice_ = som._make_lattice()
som.codebook_ = initial.copy()
som.n_features_in_ = 2
distances = np.linalg.norm(initial - x, axis=1)
c = int(np.argmin(distances))
lattice_d2 = som._lattice_squared_distances(np.array([c]))[0]
h = som._neighborhood_weights(np.array([c]), sigma)[0]
updated = initial + alpha * h[:, None] * (x - initial)
table = pd.DataFrame(
{
"unit": np.arange(1, 5),
"grid_row": som.lattice_[:, 0].astype(int),
"grid_col": som.lattice_[:, 1].astype(int),
"distance_to_x": distances,
"grid_distance_squared_to_bmu": lattice_d2,
"h_cj": h,
"m1_before": initial[:, 0],
"m2_before": initial[:, 1],
"m1_after": updated[:, 0],
"m2_after": updated[:, 1],
"is_bmu": np.arange(4) == c,
}
)
table.to_csv(output_dir / "small_example_update.csv", index=False)
(output_dir / "small_example_summary.txt").write_text(
"\n".join(
[
f"Input x = {x.tolist()}",
f"alpha = {alpha}",
f"sigma = {sigma}",
f"BMU (1-based) = {c + 1}",
"Updated codebooks:",
np.array2string(updated, precision=7),
]
),
encoding="utf-8",
)
return table
if __name__ == "__main__":
frame = run(ROOT / "outputs")
print(frame.to_string(index=False, float_format=lambda v: f"{v:.6f}"))
4.2. خروجی مورد انتظار
| واحد | مختصات شبکه | فاصله تا ورودی | فاصله شبکهای مربعشده | وزن همسایگی | بردار پیش از بهروزرسانی | بردار پس از بهروزرسانی | BMU |
| 1 | (0, 0) | 0.8246 | 1 | 0.6065 | (0.0000, 0.0000) | (0.2426, 0.0607) | خیر |
| 2 | (0, 1) | 0.2828 | 0 | 1.0000 | (1.0000, 0.0000) | (0.9000, 0.1000) | بله |
| 3 | (1, 0) | 1.1314 | 2 | 0.3679 | (0.0000, 1.0000) | (0.1472, 0.8528) | خیر |
| 4 | (1, 1) | 0.8246 | 1 | 0.6065 | (1.0000, 1.0000) | (0.9393, 0.7574) | خیر |
واحد ۲ نزدیکترین بردار کد و BMU است. وزن آن برابر یک و وزن دو همسایه مستقیم برابر است. واحد قطری با وزن حرکت میکند. خروجی عددی تا خطای گردکردن با محاسبات فصل اصلی برابر است.
.
5. مطالعه موردی اول: سازماندهی داده Wine
5.1 مسئله و داده
مجموعه Wine شامل ۱۷۸ نمونه، ۱۳ ویژگی شیمیایی و سه رقم شراب است. نسخه موجود در scikit-learn از مجموعه استاندارد UCI گرفته شده است. SOM فقط با ویژگیها آموزش میبیند. برچسب رقم شراب پس از آموزش برای دو هدف استفاده میشود:
- برچسبگذاری اکثریت واحدهای نقشه با داده آموزش و سنجش پیشبینی پسینی روی آزمون؛
- مقایسه یک خوشهبندی ثانویه بردارهای کد با برچسب واقعی.
این ارزیابی SOM را به الگوریتم نظارتشده تبدیل نمیکند، زیرا هیچ برچسبی در تابع فاصله، انتخاب BMU یا بهروزرسانی بردارهای کد وارد نشده است.
5.2 طراحی آزمایش
- تقسیم آموزش/آزمون: ۷۰/۳۰ و stratified با seed برابر ۴۱؛
- نقشه: ، شبکه مستطیلی؛
- آموزش: Batch SOM، حداکثر ۷۰ دوره؛
- مقداردهی اولیه: PCA؛
- شعاع: از ۴٫۰ تا ۰٫۸؛
- مقایسه: داده خام در برابر z-score برازششده فقط روی آموزش.
5.3 کد کامل جریان Load → Preprocess → Train → Evaluate
"""Complete Wine SOM pipeline: load -> preprocess -> train -> evaluate -> visualize."""
from __future__ import annotations
from pathlib import Path
import sys
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from sklearn.datasets import load_wine
from sklearn.metrics import (
accuracy_score,
adjusted_rand_score,
balanced_accuracy_score,
f1_score,
normalized_mutual_info_score,
)
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
ROOT = Path(__file__).resolve().parents[1]
sys.path.insert(0, str(ROOT / "src"))
from som_operational import ( # noqa: E402
SelfOrganizingMap,
majority_unit_labels,
unit_cluster_assignments,
)
def _evaluate(
som: SelfOrganizingMap,
X_train: np.ndarray,
X_test: np.ndarray,
y_train: np.ndarray,
y_test: np.ndarray,
) -> dict[str, float]:
train_bmu = som.predict_bmu(X_train)
test_bmu = som.predict_bmu(X_test)
unit_labels = majority_unit_labels(
train_bmu, y_train, som.n_units, codebook=som.codebook_
)
y_pred = unit_labels[test_bmu]
unit_clusters = unit_cluster_assignments(som.codebook_, 3, random_state=41, sample_weight=som.hit_map(X_train).ravel())
cluster_pred = unit_clusters[test_bmu]
return {
"train_qe": som.quantization_error(X_train),
"test_qe": som.quantization_error(X_test),
"train_te": som.topographic_error(X_train),
"test_te": som.topographic_error(X_test),
"unit_utilization": float(np.mean(som.hit_map(X_train).ravel() > 0)),
"posthoc_accuracy": accuracy_score(y_test, y_pred),
"posthoc_balanced_accuracy": balanced_accuracy_score(y_test, y_pred),
"posthoc_macro_f1": f1_score(y_test, y_pred, average="macro"),
"codebook_cluster_ari": adjusted_rand_score(y_test, cluster_pred),
"codebook_cluster_nmi": normalized_mutual_info_score(y_test, cluster_pred),
}
def _plot_umatrix_labels(
som: SelfOrganizingMap,
X: np.ndarray,
y: np.ndarray,
path: Path,
) -> None:
fig, ax = plt.subplots(figsize=(8, 6))
image = ax.imshow(som.u_matrix(), origin="upper", aspect="auto")
fig.colorbar(image, ax=ax, label="Mean neighbour codebook distance")
coords = som.predict_coordinates(X)
for row, col, label in zip(coords[:, 0], coords[:, 1], y, strict=True):
ax.text(col, row, str(label), ha="center", va="center", fontsize=7, alpha=0.65)
ax.set_title("Wine SOM: U-matrix with post-hoc class labels")
ax.set_xlabel("Map column")
ax.set_ylabel("Map row")
fig.tight_layout()
fig.savefig(path, dpi=180)
plt.close(fig)
def _plot_history(som: SelfOrganizingMap, path: Path) -> None:
history = som.get_history_array()
fig, ax = plt.subplots(figsize=(8, 5))
ax.plot(history[:, 0], history[:, 1], label="QE")
ax.plot(history[:, 0], history[:, 2], label="TE")
ax.set_xlabel("Epoch")
ax.set_ylabel("Metric value")
ax.set_title("Wine SOM training diagnostics")
ax.legend()
ax.grid(True, alpha=0.25)
fig.tight_layout()
fig.savefig(path, dpi=180)
plt.close(fig)
def _plot_component_planes(
som: SelfOrganizingMap, feature_names: list[str], selected: list[int], path: Path
) -> None:
planes = som.component_planes()
fig, axes = plt.subplots(1, len(selected), figsize=(4 * len(selected), 4))
axes_arr = np.atleast_1d(axes)
for ax, feature_index in zip(axes_arr, selected, strict=True):
im = ax.imshow(planes[feature_index], origin="upper", aspect="auto")
ax.set_title(feature_names[feature_index])
ax.set_xlabel("Column")
ax.set_ylabel("Row")
fig.colorbar(im, ax=ax, shrink=0.75)
fig.tight_layout()
fig.savefig(path, dpi=180)
plt.close(fig)
def run(output_dir: Path, figure_dir: Path) -> pd.DataFrame:
output_dir.mkdir(parents=True, exist_ok=True)
figure_dir.mkdir(parents=True, exist_ok=True)
data = load_wine()
X = data.data.astype(np.float64)
y = data.target.astype(np.int64)
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=0.30, stratify=y, random_state=41
)
variants: list[tuple[str, np.ndarray, np.ndarray]] = []
variants.append(("unscaled", X_train.copy(), X_test.copy()))
scaler = StandardScaler().fit(X_train)
variants.append(("zscore", scaler.transform(X_train), scaler.transform(X_test)))
rows: list[dict[str, float | str | int]] = []
trained: dict[str, tuple[SelfOrganizingMap, np.ndarray, np.ndarray]] = {}
for name, train_data, test_data in variants:
som = SelfOrganizingMap(
map_shape=(8, 6),
n_epochs=70,
training="batch",
sigma=(4.0, 0.8),
topology="rectangular",
neighborhood="gaussian",
initialization="pca",
random_state=41,
tol=1e-6,
patience=7,
batch_size=256,
record_history=True,
).fit(train_data)
metrics = _evaluate(som, train_data, test_data, y_train, y_test)
rows.append({"preprocessing": name, "n_epochs_run": som.n_iter_, **metrics})
trained[name] = (som, train_data, test_data)
result = pd.DataFrame(rows)
result.to_csv(output_dir / "wine_metrics.csv", index=False)
som, train_data, test_data = trained["zscore"]
np.save(output_dir / "wine_codebook_zscore.npy", som.codebook_)
pd.DataFrame(
som.get_history_array(),
columns=["epoch", "QE", "TE", "max_shift", "alpha", "sigma"],
).to_csv(output_dir / "wine_training_history.csv", index=False)
_plot_umatrix_labels(
som,
np.vstack([train_data, test_data]),
np.concatenate([y_train, y_test]),
figure_dir / "wine_umatrix_labels.png",
)
_plot_history(som, figure_dir / "wine_training_diagnostics.png")
_plot_component_planes(
som,
list(data.feature_names),
[0, 6, 12],
figure_dir / "wine_component_planes.png",
)
return result
if __name__ == "__main__":
metrics = run(ROOT / "outputs", ROOT / "figures")
print(metrics.to_string(index=False, float_format=lambda v: f"{v:.4f}"))
5.4 نتایج
جدول ۱ — کیفیت بازنمایی و استفاده از نقشه
| پیشپردازش | QE آموزش | QE آزمون | TE آموزش | TE آزمون | نسبت واحدهای فعال |
| بدون مقیاسبندی | 19.6239 | 19.1333 | 0.6048 | 0.6667 | 0.854 |
| استانداردسازی z-score | 1.8113 | 2.0735 | 0.0242 | 0.0741 | 0.896 |
جدول ۲ — ارزیابی پسینی با برچسبهای کنارگذاشتهشده از آموزش
| پیشپردازش | دقت متوازن | F1 کلان | ARI خوشهبندی بردارهای کد |
| بدون مقیاسبندی | 0.585 | 0.567 | 0.243 |
| استانداردسازی z-score | 0.919 | 0.925 | 0.940 |
U-matrix و برچسبهای پسینی Wine
.
منحنیهای QE و TE در آموزش Wine
.
سه صفحه مؤلفه منتخب Wine
5.5 تفسیر
مقادیر داده خام و استانداردشده مستقیماً قابل مقایسه نیستند، زیرا واحد اندازهگیری فضای ورودی عوض شده است. شواهد معنادارتر عبارتاند از:
- آزمون از حدود ۰٫۶۶۷ به ۰٫۰۷۴ کاهش یافته است؛
- دقت متوازن پسینی از حدود ۰٫۵۸۵ به ۰٫۹۱۹ رسیده است؛
- F1 کلان از ۰٫۵۶۷ به ۰٫۹۲۵ افزایش یافته است؛
- ARI خوشهبندی ثانویه بردارهای کد از حدود ۰٫۲۴۳ به ۰٫۹۴۰ رسیده است.
این اختلاف نشان میدهد که ویژگیهای با دامنه بزرگ در داده خام ساختار فاصله را غالب کردهاند. نتیجه، یک قاعده عمومی برای استانداردسازی همه دادهها نیست؛ اگر مقیاس فیزیکی نمایانگر اهمیت واقعی باشد، حفظ یا وزندهی آگاهانه آن قابل دفاع است.
5.6 عیبیابی مطالعه Wine
| نشانه | علت محتمل | اقدام پیشنهادی |
| بزرگ و نوسانی | شعاع اولیه کوچک، مقداردهی تصادفی یا مقیاس نامناسب | افزایش ، استفاده از PCA init، بررسی استانداردسازی و چند seed |
| واحدهای خالی زیاد | نقشه بیش از حد بزرگ یا داده کم | کاهش اندازه نقشه؛ گزارش نسبت واحدهای فعال |
| دقت پسینی خوب ولی U-matrix نامنظم | برچسبها در چند ناحیه پراکندهاند یا نقشه fold دارد | بررسی ، نخستین و دومین BMU و چند مقداردهی اولیه |
| آزمون بسیار بزرگتر از آموزش | تغییر توزیع، پرتهای آزمون یا نقشه کوچک | بررسی فاصلههای BMU، داده پرت و ظرفیت نقشه |
| نتیجه هر seed متفاوت | فضای داده مبهم یا آموزش ناکافی | گزارش میانگین و پراکندگی چند seed؛ Batch SOM یا PCA init |
6. مطالعه موردی دوم: نامتوازنی، نویز، داده پرت و داده ناقص
6.1 سناریو
داده مصنوعی شامل سه جمعیت با اندازههای ۵۲۰، ۱۵۰ و ۴۵ نمونه است. دو ویژگی اطلاعاتی، دو ویژگی همبسته و شش ویژگی نویزی با مقیاسهای متفاوت ساخته شدهاند. سپس ۳۵ نقطه پرت و چهار درصد مقدار مفقود افزوده شده است. برچسب 1- فقط برای ارزیابی داده پرت استفاده میشود.
سه جریان پیشپردازش مقایسه شدهاند:
- جانشانی میانه و StandardScaler؛
- جانشانی، برش صدکهای ۱ و ۹۹ و RobustScaler؛
- روش مقاوم قبلی و سپس PCA چهارمؤلفهای.
هر تنظیم با seedهای ۷، ۱۹ و ۴۳ روی یک نقشه و ۵۵ دوره Batch SOM اجرا شده است.
6.2 کد کامل مطالعه
"""Advanced SOM study with imbalance, noise, outliers, missing values, and tuning."""
from __future__ import annotations
from dataclasses import dataclass
from pathlib import Path
import sys
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from sklearn.decomposition import PCA
from sklearn.impute import SimpleImputer
from sklearn.metrics import (
adjusted_rand_score,
balanced_accuracy_score,
f1_score,
normalized_mutual_info_score,
recall_score,
roc_auc_score,
)
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import RobustScaler, StandardScaler
ROOT = Path(__file__).resolve().parents[1]
sys.path.insert(0, str(ROOT / "src"))
from som_operational import ( # noqa: E402
SelfOrganizingMap,
majority_unit_labels,
unit_cluster_assignments,
)
class QuantileClipper:
"""Train-only quantile clipping transformer with a scikit-learn-like API."""
def __init__(self, lower: float = 0.01, upper: float = 0.99) -> None:
if not 0 <= lower < upper <= 1:
raise ValueError("Quantiles must satisfy 0 <= lower < upper <= 1.")
self.lower = lower
self.upper = upper
self.lower_: np.ndarray | None = None
self.upper_: np.ndarray | None = None
def fit(self, X: np.ndarray, y: np.ndarray | None = None) -> "QuantileClipper":
arr = np.asarray(X, dtype=float)
self.lower_ = np.quantile(arr, self.lower, axis=0)
self.upper_ = np.quantile(arr, self.upper, axis=0)
return self
def transform(self, X: np.ndarray) -> np.ndarray:
if self.lower_ is None or self.upper_ is None:
raise RuntimeError("QuantileClipper is not fitted.")
return np.clip(np.asarray(X, dtype=float), self.lower_, self.upper_)
def fit_transform(self, X: np.ndarray, y: np.ndarray | None = None) -> np.ndarray:
return self.fit(X, y).transform(X)
@dataclass(frozen=True)
class SyntheticData:
X: np.ndarray
y: np.ndarray
is_outlier: np.ndarray
def make_challenging_data(seed: int = 2026) -> SyntheticData:
rng = np.random.default_rng(seed)
sizes = [520, 150, 45]
centers = [(-3.0, -1.0), (1.0, 2.6), (4.2, -2.0)]
scales = [(0.7, 0.6), (1.0, 0.8), (0.55, 1.25)]
informative = []
labels = []
for k, (n, center, scale) in enumerate(zip(sizes, centers, scales, strict=True)):
block = rng.normal(loc=center, scale=scale, size=(n, 2))
informative.append(block)
labels.append(np.full(n, k, dtype=int))
X2 = np.vstack(informative)
y = np.concatenate(labels)
# Correlated and nuisance features with deliberately different scales.
correlated = np.column_stack(
[
X2[:, 0] + rng.normal(0, 0.12, X2.shape[0]),
2.0 * X2[:, 1] + rng.normal(0, 0.25, X2.shape[0]),
]
)
noise = rng.normal(size=(X2.shape[0], 6))
noise[:, 0] *= 40.0
noise[:, 1] *= 0.02
X = np.column_stack([X2, correlated, noise])
n_outliers = 35
outliers = rng.uniform(-18.0, 18.0, size=(n_outliers, X.shape[1]))
X = np.vstack([X, outliers])
y = np.concatenate([y, np.full(n_outliers, -1, dtype=int)])
is_outlier = y == -1
# 4% missingness, introduced after outliers are added.
missing = rng.random(X.shape) < 0.04
X[missing] = np.nan
return SyntheticData(X=X, y=y, is_outlier=is_outlier)
def _split(data: SyntheticData, seed: int = 77) -> tuple[np.ndarray, ...]:
rng = np.random.default_rng(seed)
indices = np.arange(data.X.shape[0])
rng.shuffle(indices)
cut = int(0.72 * indices.size)
train_idx, test_idx = indices[:cut], indices[cut:]
return (
data.X[train_idx],
data.X[test_idx],
data.y[train_idx],
data.y[test_idx],
data.is_outlier[train_idx],
data.is_outlier[test_idx],
)
def _preprocessors() -> dict[str, Pipeline]:
return {
"standard": Pipeline(
[
("imputer", SimpleImputer(strategy="median")),
("scaler", StandardScaler()),
]
),
"robust_clip": Pipeline(
[
("imputer", SimpleImputer(strategy="median")),
("clip", QuantileClipper(0.01, 0.99)),
("scaler", RobustScaler(quantile_range=(10.0, 90.0))),
]
),
"robust_clip_pca": Pipeline(
[
("imputer", SimpleImputer(strategy="median")),
("clip", QuantileClipper(0.01, 0.99)),
("scaler", RobustScaler(quantile_range=(10.0, 90.0))),
("pca", PCA(n_components=4, random_state=77)),
]
),
}
def _evaluate_run(
som: SelfOrganizingMap,
X_train: np.ndarray,
X_test: np.ndarray,
y_train: np.ndarray,
y_test: np.ndarray,
outlier_train: np.ndarray,
outlier_test: np.ndarray,
) -> dict[str, float]:
# Labels are excluded from SOM fitting. Outliers are excluded only from the
# majority-label table so they cannot become a fictitious fourth class.
normal_train = ~outlier_train
normal_test = ~outlier_test
train_bmu = som.predict_bmu(X_train)
test_bmu = som.predict_bmu(X_test)
unit_labels = majority_unit_labels(
train_bmu[normal_train],
y_train[normal_train],
som.n_units,
codebook=som.codebook_,
)
y_pred = unit_labels[test_bmu[normal_test]]
unit_clusters = unit_cluster_assignments(som.codebook_, 3, random_state=77, sample_weight=som.hit_map(X_train).ravel())
cluster_pred = unit_clusters[test_bmu[normal_test]]
qdist = som.quantization_distances(X_test)
minority_label = 2
return {
"train_qe": som.quantization_error(X_train),
"test_qe": som.quantization_error(X_test),
"train_te": som.topographic_error(X_train),
"test_te": som.topographic_error(X_test),
"unit_utilization": float(np.mean(som.hit_map(X_train).ravel() > 0)),
"posthoc_balanced_accuracy": balanced_accuracy_score(y_test[normal_test], y_pred),
"posthoc_macro_f1": f1_score(y_test[normal_test], y_pred, average="macro"),
"minority_recall": recall_score(
y_test[normal_test], y_pred, labels=[minority_label], average="macro", zero_division=0
),
"codebook_cluster_ari": adjusted_rand_score(y_test[normal_test], cluster_pred),
"codebook_cluster_nmi": normalized_mutual_info_score(y_test[normal_test], cluster_pred),
"outlier_auc_from_bmu_distance": roc_auc_score(outlier_test.astype(int), qdist),
}
def _plot_summary(summary: pd.DataFrame, path: Path) -> None:
metrics = ["test_qe_mean", "test_te_mean", "posthoc_balanced_accuracy_mean", "outlier_auc_from_bmu_distance_mean"]
labels = ["Test QE", "Test TE", "Balanced accuracy", "Outlier AUC"]
fig, axes = plt.subplots(2, 2, figsize=(10, 8))
for ax, metric, label in zip(axes.ravel(), metrics, labels, strict=True):
ax.bar(summary["preprocessing"], summary[metric])
ax.set_title(label)
ax.tick_params(axis="x", rotation=25)
ax.grid(True, axis="y", alpha=0.2)
fig.tight_layout()
fig.savefig(path, dpi=180)
plt.close(fig)
def _plot_selected_map(
som: SelfOrganizingMap,
X_test: np.ndarray,
y_test: np.ndarray,
outlier_test: np.ndarray,
path: Path,
) -> None:
fig, ax = plt.subplots(figsize=(8, 6))
image = ax.imshow(som.u_matrix(), origin="upper", aspect="auto")
fig.colorbar(image, ax=ax, label="U-matrix")
coords = som.predict_coordinates(X_test)
for (row, col), label, is_outlier in zip(coords, y_test, outlier_test, strict=True):
text = "O" if is_outlier else str(label)
ax.text(col, row, text, ha="center", va="center", fontsize=6, alpha=0.55)
ax.set_title("Robust preprocessing: classes and outliers on the SOM")
ax.set_xlabel("Map column")
ax.set_ylabel("Map row")
fig.tight_layout()
fig.savefig(path, dpi=180)
plt.close(fig)
def run(output_dir: Path, figure_dir: Path) -> tuple[pd.DataFrame, pd.DataFrame]:
output_dir.mkdir(parents=True, exist_ok=True)
figure_dir.mkdir(parents=True, exist_ok=True)
data = make_challenging_data()
X_train_raw, X_test_raw, y_train, y_test, out_train, out_test = _split(data)
raw_rows: list[dict[str, float | int | str]] = []
fitted_for_plot: tuple[SelfOrganizingMap, np.ndarray] | None = None
for name, pipeline in _preprocessors().items():
X_train = pipeline.fit_transform(X_train_raw)
X_test = pipeline.transform(X_test_raw)
for seed in (7, 19, 43):
som = SelfOrganizingMap(
map_shape=(7, 5),
n_epochs=55,
training="batch",
sigma=(3.5, 0.8),
topology="rectangular",
neighborhood="gaussian",
initialization="random_samples",
random_state=seed,
tol=1e-5,
patience=6,
batch_size=256,
record_history=False,
).fit(X_train)
metrics = _evaluate_run(
som, X_train, X_test, y_train, y_test, out_train, out_test
)
raw_rows.append(
{
"preprocessing": name,
"seed": seed,
"n_features_after_preprocessing": X_train.shape[1],
"n_epochs_run": som.n_iter_,
**metrics,
}
)
if name == "robust_clip_pca" and seed == 19:
fitted_for_plot = (som, X_test)
raw = pd.DataFrame(raw_rows)
raw.to_csv(output_dir / "robustness_runs.csv", index=False)
metric_columns = [
c for c in raw.columns if c not in {"preprocessing", "seed", "n_features_after_preprocessing"}
]
summary_parts = []
for name, group in raw.groupby("preprocessing", sort=False):
row: dict[str, float | str] = {"preprocessing": name}
for metric in metric_columns:
row[f"{metric}_mean"] = float(group[metric].mean())
row[f"{metric}_std"] = float(group[metric].std(ddof=1))
summary_parts.append(row)
summary = pd.DataFrame(summary_parts)
summary.to_csv(output_dir / "robustness_summary.csv", index=False)
_plot_summary(summary, figure_dir / "robustness_comparison.png")
if fitted_for_plot is not None:
som, X_test = fitted_for_plot
_plot_selected_map(
som, X_test, y_test, out_test, figure_dir / "robustness_umatrix.png"
)
return raw, summary
if __name__ == "__main__":
_, summary_frame = run(ROOT / "outputs", ROOT / "figures")
print(summary_frame.to_string(index=False, float_format=lambda v: f"{v:.4f}"))
6.3 نتایج چند seed
جدول ۳ — کیفیت نقشه در چند seed؛ میانگین ± انحراف معیار
| پیشپردازش | QE آزمون | TE آزمون | نسبت واحدهای فعال |
| جانشانی + StandardScaler | 1.4281 ± 0.0058 | 0.1000 ± 0.0082 | 0.990 |
| جانشانی + برش صدکی + RobustScaler | 1.5160 ± 0.0444 | 0.1175 ± 0.0317 | 0.800 |
| روش مقاوم + PCA چهار مؤلفهای | 1.1125 ± 0.0351 | 0.1460 ± 0.0458 | 0.800 |
جدول ۴ — ارزیابی پسینی و تشخیص داده پرت
| پیشپردازش | دقت متوازن | F1 کلان | بازخوانی اقلیت | AUC فاصله BMU |
| جانشانی + StandardScaler | 0.927 ± 0.032 | 0.946 ± 0.025 | 0.796 ± 0.085 | 1.000 |
| جانشانی + برش صدکی + RobustScaler | 0.935 ± 0.002 | 0.954 ± 0.002 | 0.852 ± 0.032 | 1.000 |
| روش مقاوم + PCA چهار مؤلفهای | 0.652 ± 0.017 | 0.611 ± 0.017 | 0.019 ± 0.032 | 0.996 |
مقایسه پیشپردازشها در مطالعه مقاومت
.
U-matrix نمونه روش مقاوم و PCA
.
6.4 تحلیل تفاوتها
برش صدکی و مقیاسبندی مقاوم دقت متوازن را از ۰٫۹۲۷ به ۰٫۹۳۵ و بازخوانی کلاس اقلیت را از ۰٫۷۹۶ به ۰٫۸۵۲ افزایش دادهاند. این بهبود با افزایش و کاهش نسبت واحدهای فعال از حدود ۰٫۹۹ به ۰٫۸۰ همراه است. بنابراین پیشپردازش مقاوم، یک بهبود بدون هزینه نیست؛ بخشهایی از نقشه فشردهتر و کماستفادهتر شدهاند.
PCA چهارمؤلفهای عددی کوچکتری ایجاد کرده است، اما این عدد با حالتهای دهویژگی قابل مقایسه نیست. مهمتر آنکه بازخوانی کلاس اقلیت به حدود ۰٫۰۱۹ سقوط کرده است. علت محتمل آن است که مؤلفههای بیشینهواریانس، ساختار جمعیت کوچک را کماهمیت کردهاند. این نتیجه نشان میدهد کاهش بُعد بدون نظارت ممکن است ساختار نادر را حذف کند، حتی اگر خطای کوانتیزهسازی در فضای جدید مطلوب به نظر برسد.
AUC فاصله BMU برای تشخیص نقاط پرت نزدیک یک است، زیرا نقاط پرت این آزمایش عمداً فاصله زیادی از جمعیت اصلی دارند. این نتیجه یک بررسی سلامت pipeline است و نباید بهعنوان معیار عملکرد عمومی تشخیص ناهنجاری SOM گزارش شود.
در دو پیشپردازش مقاوم، K-means ثانویه روی بردارهای کد نتوانست افراز برچسبی قابلقبولی بسازد، در حالی که برچسبگذاری محلی واحدها عملکرد خوبی داشت. این اختلاف یادآوری میکند که «واحد نقشه»، «ناحیه توپولوژیک» و «خوشه K-means روی codebook» سه شیء متفاوتاند.
.
6.5 عیبیابی پیشرفته و تنظیم
| مشکل | تشخیص کمی | مداخله | خطر مداخله |
| کشیدهشدن prototypeها به سمت پرتها | صدکهای بالای فاصله BMU و codebookهای دور | RobustScaler، برش دامنهای یا مدل مقاوم | حذف یا تضعیف رخداد واقعی نادر |
| از دست رفتن کلاس اقلیت | بازخوانی اقلیت پایین با دقت کلی مناسب | افزایش ظرفیت، نمونهبرداری کنترلشده فقط در جریان آموزش آنلاین، حفظ ویژگیهای نادر | تغییر توزیع مؤثر آموزش |
| زیاد | نخستین و دومین BMU اغلب غیرمجاور | افزایش شعاع اولیه، کاهش آهستهتر، PCA init، تغییر نسبت نقشه | افزایش بیشهموارسازی و |
| واحد مرده | hit map صفر برای بخش بزرگی از نقشه | کوچککردن نقشه یا افزایش شعاع آغازین | کاهش تفکیک محلی |
| PCA نتیجه را بدتر میکند | افت اقلیت یا ساختار برچسبی با وجود کمتر | افزایش مؤلفهها، تحلیل بارها، حذف PCA | بازگشت هزینه ابعادی و نویز |
| ناپایداری seed | انحراف معیار زیاد QE/TE یا metric پسینی | اجرای حداقل ۱۰ seed در پژوهش نهایی | هزینه محاسباتی بیشتر |
| حافظه زیاد در Batch SOM | رشد در پیادهسازی ساده | محاسبه بلوکی BMU و انباشت صورت/مخرج | افزایش تعداد عملیات و I/O |
| داده ناقص جدید در زمان پیشبینی | خطای صریح ورودی | ذخیره و استفاده از imputer آموزش | رانش الگوی فقدان باید جداگانه پایش شود |
6.6 راهبرد تنظیم پیشنهادی
- یک pipeline پایه با جانشانی و مقیاسبندی بسازید.
- چند اندازه نقشه را با تعداد واحدهای ثابتنسبت به داده بررسی کنید.
- برای هر اندازه، دستکم چند مقدار و را مقایسه کنید.
- همه تنظیمها را با seedهای متعدد اجرا کنید.
- ، ، نسبت واحدهای فعال و پایداری U-matrix را همزمان گزارش کنید.
- اگر برچسب مستقل وجود دارد، فقط پس از آموزش برای ارزیابی استفاده کنید.
- پیشپردازشهایی که فضای ویژگی متفاوت میسازند را صرفاً بر اساس مقدار خام رتبهبندی نکنید.
- سادهترین نقشهای را انتخاب کنید که نتیجه پایدار و تفسیرپذیر ارائه دهد.
7. جمعبندی نهایی
7.1 نکات اجرایی اصلی
- SOM پایه باید ورودی عددی کامل و مقیاسبندیشده بر اساس منطق دامنه دریافت کند.
- نسخه آنلاین به ترتیب نمونه و seed حساس است؛ نسخه دستهای نسبت به جایگشت سطرها ناورداست، ولی همچنان به مقداردهی اولیه و تساویشکنی وابسته است.
- تعریف ، و مجاورت شبکه باید همراه هر نتیجه ثبت شود.
- اندازه نقشه، شعاع اولیه و روش مقداردهی اولیه سه عامل اصلی نظم توپولوژیکاند.
- U-matrix و hit map ابزار تشخیصاند؛ تعداد نواحی بصری آنها تعداد خوشه قطعی نیست.
- فاصله BMU میتواند برای امتیاز ناهنجاری مفید باشد، اما آستانه باید روی داده اعتبارسنجی تعیین شود.
- کاهش بُعد یا معیار فاصله جدید یک جزئیات خنثی نیست و میتواند ساختار مسئله را عوض کند.
7.2 هشدارهای پیادهسازی
- NaN را داخل تابع فاصله نادیده نگیرید، مگر اینکه نرمالسازی بر تعداد ویژگیهای مشترک و پیامدهای آن را صریح تعریف کرده باشید.
- ماتریس کامل فاصله را برای داده بزرگ نگه ندارید؛ BMU را بلوکی محاسبه کنید.
- در آموزش آنلاین، کاهش بسیار سریع میتواند fold را تثبیت کند.
- QE کمتر روی نقشه بزرگتر یا فضای PCAشده الزاماً مدل بهتر را نشان نمیدهد.
- seed، نسخه کتابخانهها و قرارداد تساوی BMU باید در مخزن نتیجه ثبت شوند.
- برچسبگذاری اکثریت واحدها یک لایه ارزیابی یا کاربرد پسینی است و نباید بهعنوان بخش SOM بدون نظارت معرفی شود.
7.3 توسعههای عملیاتی بعدی
- پیادهسازی شبکه ششضلعی و تعریف متناظر و U-matrix؛
- افزودن partial_fit همراه با فراموشی برای جریان دارای رانش؛
- جستوجوی تقریبی BMU برای نقشههای بسیار بزرگ؛
- SOM مقاوم با prototypeهای medoid یا تابع زیان heavy-tailed؛
- نسخه GPU با PyTorch، همراه با آزمون برابری عددی؛
- prototype و بهروزرسانی سازگار با DTW برای سری زمانی؛
- ارزیابی پایداری نقشه با bootstrap و تطبیق واحدهای اجراهای مستقل.
7.4 اقلام بسته بازتولیدپذیری
- src/som_operational/som.py: پیادهسازی کامل؛
- case_studies/small_example.py: مثال دقیق فصل؛
- case_studies/wine_case_study.py: مطالعه Wine؛
- case_studies/robustness_case_study.py: مطالعه نویز و عدمتوازن؛
- tests/test_som.py: ده آزمون خودکار؛
- outputs/*.csv: جدولهای خام؛
- figures/*.png: نمودارهای گزارش؛
- outputs/environment.json: نسخهها و seedها؛
- outputs/sha256.json: checksum خروجیها؛
- requirements-tested.txt: محیط دقیق؛
- requirements-current-stable.txt: نسخههای پایدار بررسیشده.
منابع فنی و داده
Kohonen, T. (2001). Self-Organizing Maps (3rd ed.). Springer.
UCI Machine Learning Repository. Wine Data Set.
scikit-learn developers. The Wine recognition dataset and preprocessing documentation.
NumPy, SciPy, pandas, Matplotlib, pytest, Python and scikit-learn official release documentation, accessed 2026-08-05.
