Self-Organizing Maps: Complete Guide — Kohonen Maps and Topology

Summary

Self-Organizing Maps (SOM), or self-organizing maps, are one of the most elegant unsupervised learning algorithms ever designed. Introduced by Finnish professor Teuvo Kohonen in 1982, these artificial neural networks have a remarkable property: they project high-dimensional data onto a two-dimensional map while preserving the topology of the input space. In other words, similar samples in the original space end up close to each other on the map. This characteristic makes SOM a powerful tool for multidimensional data visualization, clustering, and exploring hidden structures in complex datasets. Unlike traditional clustering methods such as K-Means, a SOM doesn’t simply produce isolated groups: it builds an organized spatial representation where geographical proximity on the map reflects similarity between the data.

Mathematical Principle of the Algorithm

The operation of a Self-Organizing Map is based on a competitive learning mechanism structured around three fundamental steps. Each step is essential to understanding how the map manages to organize itself.

Step 0: Initialization

A grid of neurons (usually rectangular or hexagonal) is set up. Each neuron i has a weight vector wᵢ of the same dimension as the input data. These weights are initialized randomly, meaning that at the start, the map is completely disorganized and does not reflect any structure of the data.

Step 1: Finding the Best Matching Unit (BMU)

For each training sample x, the neuron whose weight is closest to x is found according to a distance metric, typically the Euclidean distance:

$$c = \arg\min_i ||\mathbf{x} – \mathbf{w}_i||$$

The neuron c thus identified is called the Best Matching Unit (BMU), i.e., the best-matching unit. It is the neuron whose weight vector most closely resembles the current sample. This competition between neurons is at the heart of SOM learning.

Step 2: Weight Update

Once the BMU is identified, not only its weights are updated, but also those of its neighbors on the grid. The update rule is as follows:

$$\mathbf{w}_i(t+1) = \mathbf{w}_i(t) + \alpha(t) \cdot h(c, i, t) \cdot (\mathbf{x} – \mathbf{w}_i(t))$$

where:

• α(t) is the learning rate at iteration t, which decreases over time.
• h(c, i, t) is the Gaussian neighborhood function that measures the influence of the BMU c on neighboring neuron i.

Step 3: Neighborhood Function

The neighborhood function determines how a neuron i is affected by the presence of the BMU c:

$$h(c, i, t) = \exp\left(-\frac{||\mathbf{r}_c – \mathbf{r}_i||^2}{2 \cdot \sigma(t)^2}\right)$$

where r_c and r_i are the respective positions of neurons c and i on the two-dimensional grid. The Gaussian function ensures that neurons close to the BMU are strongly updated, while distant neurons are updated very little.

Parameter Decay

To ensure convergence of the algorithm, two parameters decay exponentially over the iterations:

$$\alpha(t) = \alpha_0 \cdot \exp(-t / \text{max_iterations})$$

$$\sigma(t) = \sigma_0 \cdot \exp(-t / \text{max_iterations})$$

• α₀ is the initial learning rate.
• σ₀ is the initial neighborhood radius.

This decay is crucial: in the early iterations, a large neighborhood radius allows the map to organize globally. Progressively, the radius shrinks and learning becomes more local, refining the details of the map. This is the process that allows the SOM to preserve topology: neurons that are close on the grid end up with similar weights, thus reflecting the structure of the input data.

Intuition: The Archipelago of Ecosystems

Imagine an archipelago of islands arranged on a regular grid in the middle of the ocean. Each island represents a neuron of the self-organizing map. At the beginning, each island has a completely random habitat: some are desert, some are covered in tropical forests, some are snowy, others are swamps. No order is visible from the sky.

Now imagine that animals representing your data arrive progressively, one by one. Each animal seeks the island whose habitat matches it best. A polar bear will naturally choose the snowy island, a toucan will prefer the tropical island. This elected island is the Best Matching Unit (BMU) for that animal.

But here is the magic of the algorithm: when an island receives a visitor, it adapts its habitat to even better match the animal, and its neighboring islands do the same, but to a lesser extent. The island itself changes a lot, its direct neighbors change a little, and the more distant islands barely change at all. This is the neighborhood function in action.

Over time, with thousands of animals visiting the archipelago, something remarkable happens: islands close to each other end up with similar habitats. The archipelago organizes itself spontaneously! Snowy islands cluster together in one corner, tropical forests in another, deserts in a third. The map becomes an organized representation of the ecosystems present in the data. Similar animals always end up in the same region of the archipelago. Topology is preserved: if you walk across the map, you pass through ecosystems that change gradually, without any abrupt jumps.

This is exactly what a Self-Organizing Map does with your data: it builds a map where spatial proximity reflects semantic similarity.

Python Implementation

From-Scratch Implementation of a 2D SOM

Here is a complete implementation of a Self-Organizing Map from scratch, with no dependencies other than NumPy and Matplotlib. This version uses a rectangular grid and a Gaussian neighborhood function.

import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import load_iris
from sklearn.preprocessing import MinMaxScaler

class SelfOrganizingMap:
    """
    From-scratch implementation of a Kohonen Self-Organizing Map (SOM)
    with a 2D rectangular grid.
    """

    def __init__(self, map_width, map_height, input_dim,
                 sigma0=None, learning_rate0=0.5, n_iterations=10000,
                 seed=42):
        self.map_width = map_width
        self.map_height = map_height
        self.input_dim = input_dim
        self.n_iterations = n_iterations
        self.learning_rate0 = learning_rate0
        # Initial radius: max of half the map dimensions
        self.sigma0 = sigma0 if sigma0 else max(map_width, map_height) / 2.0

        rng = np.random.RandomState(seed)
        # Random weight initialization
        self.weights = rng.rand(map_width, map_height, input_dim)

        # Position grid (2D coordinates)
        self.grid_x = np.arange(map_width)
        self.grid_y = np.arange(map_height)
        self.pos_x, self.pos_y = np.meshgrid(self.grid_x, self.grid_y)
        self.positions = np.column_stack([self.pos_x.ravel(), self.pos_y.ravel()])

    def _compute_distances(self, x):
        """Computes the Euclidean distance between x and all weights."""
        diff = self.weights - x[np.newaxis, np.newaxis, :]
        return np.sqrt(np.sum(diff ** 2, axis=2))

    def find_bmu(self, x):
        """Finds the Best Matching Unit for a sample x."""
        distances = self._compute_distances(x)
        idx = np.argmin(distances)
        bmu_x, bmu_y = np.unravel_index(idx, (self.map_width, self.map_height))
        return bmu_x, bmu_y

    def _gaussian_neighborhood(self, bmu_x, bmu_y, sigma):
        """Computes the Gaussian neighborhood function for all neurons."""
        dist_sq = ((self.pos_x - bmu_x) ** 2 + (self.pos_y - bmu_y) ** 2)
        neighborhood = np.exp(-dist_sq / (2 * sigma ** 2))
        return neighborhood

    def fit(self, X, verbose=True):
        """Trains the SOM on data X."""
        for t in range(self.n_iterations):
            # Random sample selection
            idx = np.random.randint(0, len(X))
            x = X[idx]

            # Parameter decay
            alpha = self.learning_rate0 * np.exp(-t / self.n_iterations)
            sigma = self.sigma0 * np.exp(-t / self.n_iterations)

            # Find the BMU
            bmu_x, bmu_y = self.find_bmu(x)

            # Neighborhood function
            neighborhood = self._gaussian_neighborhood(bmu_x, bmu_y, sigma)
            neighborhood = neighborhood[:, :, np.newaxis]  # Broadcasting

            # Weight update
            self.weights += alpha * neighborhood * (x - self.weights)

            if verbose and t % 1000 == 0:
                print(f"Iteration {t}/{self.n_iterations} | "
                      f"α={alpha:.4f} | σ={sigma:.4f}")

        return self

    def transform(self, X):
        """Projects data X onto the map (returns BMU coordinates)."""
        bmu_coords = []
        for x in X:
            bx, by = self.find_bmu(x)
            bmu_coords.append((bx, by))
        return np.array(bmu_coords)

    def u_matrix(self):
        """Computes the U-Matrix (map of mean distances between neighbors)."""
        u = np.zeros((self.map_width, self.map_height))
        for i in range(self.map_width):
            for j in range(self.map_height):
                neighbors = []
                for di, dj in [(-1,0),(1,0),(0,-1),(0,1)]:
                    ni, nj = i+di, j+dj
                    if 0 <= ni < self.map_width and 0 <= nj < self.map_height:
                        dist = np.linalg.norm(self.weights[i,j] - self.weights[ni,nj])
                        neighbors.append(dist)
                if neighbors:
                    u[i,j] = np.mean(neighbors)
                else:
                    u[i,j] = 0
        return u

Application on the Iris Dataset

Let’s now illustrate the use of our SOM on Fisher’s famous Iris dataset, a machine learning classic containing 150 samples of three iris species.

# Load and normalize data
iris = load_iris()
X = iris.data
y = iris.target
scaler = MinMaxScaler()
X_scaled = scaler.fit_transform(X)

# Create and train the SOM
som = SelfOrganizingMap(
    map_width=10,
    map_height=10,
    input_dim=4,
    sigma0=4.0,
    learning_rate0=0.5,
    n_iterations=50000,
    seed=42
)
som.fit(X_scaled, verbose=True)

# U-Matrix visualization
u_mat = som.u_matrix()
plt.figure(figsize=(8, 6))
plt.imshow(u_mat.T, cmap='viridis')
plt.colorbar(label="Mean distance between neighboring neurons")
plt.title("U-Matrix of the SOM applied to the Iris dataset")
plt.xlabel("X coordinate on the grid")
plt.ylabel("Y coordinate on the grid")
plt.tight_layout()
plt.savefig("som_iris_umatrix.png", dpi=150)
plt.show()

# Mark the BMU of each class on the grid
bmu_coords = som.transform(X_scaled)
plt.figure(figsize=(8, 6))
colors = ['red', 'green', 'blue']
species = ['Setosa', 'Versicolor', 'Virginica']
for cls in range(3):
    mask = y == cls
    bx, by = bmu_coords[mask].T
    plt.scatter(bx, by, c=colors[cls], label=species[cls], alpha=0.6,
                edgecolors='black', s=100)
plt.title("SOM Projection — Iris Species on the Grid")
plt.xlabel("X coordinate")
plt.ylabel("Y coordinate")
plt.legend()
plt.gca().set_xlim(-1, som.map_width)
plt.gca().set_ylim(-1, som.map_height)
plt.tight_layout()
plt.savefig("som_iris_clusters.png", dpi=150)
plt.show()

Using MiniSom

For practical applications, the MiniSom library provides a robust and optimized implementation of SOMs. It can be installed easily via pip.

pip install minisom

from minisom import MiniSom

# MiniSom configuration
som_minisom = MiniSom(
    x=10, y=10,
    input_len=4,
    sigma=4.0,
    learning_rate=0.5,
    neighborhood_function='gaussian'
)

# Random initialization
som_minisom.random_weights_init(X_scaled)

# Training with decay
som_minisom.train_random(
    X_scaled,
    num_iteration=50000,
    verbose=True
)

# Weight map visualization
plt.figure(figsize=(10, 8))
for i in range(4):
    plt.subplot(2, 2, i+1)
    plt.imshow(som_minisom.get_weights()[:, :, i].T, cmap='plasma')
    feature_names = ['Sepal length', 'Sepal width', 'Petal length', 'Petal width']
    plt.title(f"Weight Map — {feature_names[i]}")
    plt.colorbar()
plt.tight_layout()
plt.savefig("som_weight_maps.png", dpi=150)
plt.show()

RGB Color Clustering with SOM

A particularly visual example of SOM usefulness is color clustering: take an image, extract all its RGB pixels, and train a SOM to group them into a representative palette.

from skimage import io

# Load an image (replace with your path)
image = io.imread("photo.jpg")
pixels = image.reshape(-1, 3).astype(np.float64) / 255.0

# SOM to extract a palette of 16 colors (4x4)
color_som = SelfOrganizingMap(
    map_width=4,
    map_height=4,
    input_dim=3,
    sigma0=2.0,
    learning_rate0=0.5,
    n_iterations=20000,
    seed=123
)
color_som.fit(pixels, verbose=False)

palette = color_som.weights.reshape(4, 4, 3)

# Display the palette
fig, axes = plt.subplots(1, 2, figsize=(12, 5))

axes[0].imshow(image)
axes[0].set_title("Original image")
axes[0].axis("off")

axes[1].imshow(palette)
axes[1].set_title("SOM Palette (4×4)")
axes[1].axis("off")

plt.tight_layout()
plt.savefig("som_color_palette.png", dpi=150)
plt.show()

# Reconstruct the image with the reduced palette
from scipy.spatial.distance import cdist
distances = cdist(pixels, color_som.weights.reshape(16, 3))
assigned = np.argmin(distances, axis=1)
reconstructed = color_som.weights.reshape(16, 3)[assigned].reshape(image.shape)

plt.imshow(reconstructed)
plt.title(f"Image reconstructed with {16} SOM colors")
plt.axis("off")
plt.savefig("som_reconstructed.png", dpi=150)
plt.show()

This application perfectly demonstrates the power of SOMs: the map has automatically extracted a palette of 16 dominant colors that preserves the chromatic distribution of the original image, while organizing these colors in a topological manner on the grid.

Key Hyperparameters

The choice of hyperparameters is critical for the quality of a Self-Organizing Map. Here are the main ones:

Practical Tips

• Map size: a rule of thumb suggests approximately 5 × √N neurons, where N is the number of samples. For 1,000 data points, a 15×15 grid (225 neurons) is appropriate.
• σ₀ too small: the map fragments and doesn’t organize globally.
• σ₀ too large: organization is too slow, and local details are poorly captured.
• Learning rate too high: instability and weight oscillation.
• Learning rate too low: very slow convergence, suboptimal map.
• Number of iterations: sufficient iterations are needed for σ(t) and α(t) to decay to very low values, ensuring final stabilization.
• Data normalization: always normalize data before training (MinMaxScaler or StandardScaler). SOMs are sensitive to different scales between variables.

Advantages and Limitations

Advantages

1. Topology preservation: the distinguishing feature of SOMs. Similar data points are close on the map, providing an intuitive visualization of high-dimensional structures.
2. Intuitive visualization: the U-Matrix clearly reveals the boundaries between clusters, where distances between neighboring neurons are high.
3. Unsupervised and assumption-free: unlike K-Means, no assumptions about cluster shape are necessary (sphericity, equal size, etc.).
4. Robustness to noise: the collective nature of learning (neighborhood update) makes SOMs less sensitive to outliers than individual methods.
5. Natural cluster detection: the map reveals the underlying data structure without requiring the number of clusters to be specified in advance.
6. Interpretability: the weight map can be visualized for each individual feature, providing a fine-grained understanding of what each region of the map represents.

Limitations

1. Non-determinism: the random weight initialization and sampling order make results unreproducible without setting a random seed.
2. Tricky hyperparameter selection: map size, sigma, learning rate, and number of iterations often require experience and multiple trials.
3. Computational cost: each iteration requires computing the distance between the sample and all weights. For very large maps and lots of data, training can be long.
4. No intrinsic quality metric: unlike K-Means which minimizes intra-cluster inertia, SOMs have no single objective score to evaluate result quality.
5. Difficult to evaluate quantitatively: visualization is subjective, and there are few standard metrics to compare two maps.
6. Normalization sensitivity: non-normalized variables dominate the Euclidean distance and completely bias the result.

4 Concrete Use Cases

1. Customer Data Exploration in Marketing

SOMs are excellent for segmenting a customer base without prior assumptions about the number of segments. By projecting demographic, behavioral, and transactional data onto a 2D map, natural zones emerge: premium customers, occasional customers, new users, etc. The topology of the map also reveals transitions between segments, identifying customers on the boundary who could shift from one behavior to another.

2. Document and Text Analysis

By representing each document as a TF-IDF vector or embeddings (Word2Vec, BERT), a SOM can automatically organize a document collection. Articles on similar topics cluster together on the map. A journalist could thus explore archives of thousands of articles by visually navigating the map, identifying emerging themes and outlier articles that deserve attention.

3. Medical Image Diagnostic

SOMs can analyze features extracted from medical images (X-rays, MRIs, CT scans) to detect anomalies. By training a map on healthy tissue, regions corresponding to pathological tissue appear as distinct areas of the U-Matrix, facilitating visual identification of suspicious areas by radiologists.

4. Industrial Quality Control

In manufacturing, SOMs monitor hundreds of sensors simultaneously. Normal operation occupies a specific region of the map. The appearance of samples in new areas signals a drift or emerging defect. This approach detects multivariate anomalies that univariate thresholds on each individual sensor would completely miss.

See Also

• Is it true that I can buy an iPhone with only 1 euro?
• Master Grid Point Calculation in Lattice Cubes with Python: Complete Guide and Tutorials