Image Segmentation: from Pixel Clustering to Encoder–Decoders

Segmentation assigns a label to every pixel.

The lectures distinguish problems by what that label means:

• Unlabelled regions. The label is only an ID: pixels with the same ID belong together. This is what clustering, superpixels and graph methods produce (lecture 19).
• Semantic segmentation. The label is a category such as cow, grass or sky. Two cows get the same label.
• Instance segmentation. Each countable object gets its own ID: cow 1, cow 2.
• Panoptic segmentation. Both at once: every pixel gets a category, and pixels of countable objects also get an instance ID.¹

Slide 41 adds the vocabulary for this split. Things are countable objects (cow, person, car). Stuff is regions or materials (road, grass, sky).² Instance segmentation focuses on things.

The simplest segmenter ignores where pixels are and looks only at their colour.

Each pixel becomes a point in colour space, for example \((R, G, B)\) or \((L, a, b)\). Clusters of points become segments.

k-means looks for \(k\) centres \(\boldsymbol\mu_1, \dots, \boldsymbol\mu_k\) and an assignment \(c_i\) of each point that minimise the within-cluster sum of squares:³

$$J = \sum_{i=1}^{n} \big\lVert \mathbf x_i - \boldsymbol\mu_{c_i} \big\rVert^2 .$$

Lloyd's algorithm alternates two steps until nothing changes:

1. Assignment step. Give each point the label of its nearest centre.
2. Update step. Move each centre to the mean of the points assigned to it.

Slide 12 replaces hard assignments with a Gaussian mixture model. Each cluster becomes a Gaussian with its own covariance, and each pixel gets a probability of belonging to each cluster. k-means is the limiting case where every cluster has the same tiny spherical covariance.

Mean shift treats the points as samples from an unknown probability density, and looks for its peaks (modes).⁴

Each point is then labelled by the peak it climbs to. The number of clusters is not chosen in advance: it is the number of peaks.

The procedure on slide 24 is:

1. Centre a window of radius \(h\) (the bandwidth) on a point.
2. Compute the mean of the data inside the window.
3. Move the window to that mean, and repeat until it stops moving.
4. Points whose windows end at the same peak form one cluster.

The move from \(\mathbf x\) to the weighted mean is the mean-shift vector (slide 23). The weights come from a kernel profile \(g\):

$$\mathbf m(\mathbf x) = \frac{\displaystyle\sum_{i=1}^{n} \mathbf x_i\, g\!\left(\left\lVert \frac{\mathbf x - \mathbf x_i}{h}\right\rVert^2\right)}{\displaystyle\sum_{i=1}^{n} g\!\left(\left\lVert \frac{\mathbf x - \mathbf x_i}{h}\right\rVert^2\right)} - \mathbf x .$$

A flat kernel gives \(g = 1\) inside the window and \(0\) outside, so the update is the plain mean of the points in the window. A Gaussian kernel weights every point by its distance, with \(g(u) = e^{-u/2}\).

Clustering in colour and position. Slide 25 appends the pixel coordinates to the colour, giving five-dimensional features \((x, y, L, a, b)\).

Each dimension has its own bandwidth, or equivalently, the coordinates are rescaled before one shared bandwidth is used (the notebook does this in cell 25). A larger spatial bandwidth relative to the colour bandwidth gives smoother, more blob-like regions (slide 26).

Over-segmentation splits the image into regions that are smaller than objects. Objects are still separated from the background, but each object is broken into many parts.

These parts are superpixels: groups of adjacent pixels with similar colour. They make a compact intermediate representation: a \(600 \times 800\) image of \(480\,000\) pixels becomes about \(200\) superpixels (slide 35).

SLIC (simple linear iterative clustering) is k-means in \((l, a, b, x, y)\) space with two changes:⁵

1. Initialisation. Place \(k\) centres on a regular grid with spacing \(S\), then move each to the lowest-gradient pixel in its \(3 \times 3\) neighbourhood, so it does not start on an edge.
2. Local search. Each centre only competes for pixels in a \(2S \times 2S\) window around it. This makes the cost linear in the number of pixels.

The distance mixes colour and position (slide 32):

$$D = d_{lab} + \frac{\alpha}{S}\, d_{xy},$$ $$d_{lab} = \sqrt{(l - l_k)^2 + (a - a_k)^2 + (b - b_k)^2}, \qquad d_{xy} = \sqrt{(x - x_k)^2 + (y - y_k)^2}.$$

The weight \(\alpha\) is the compactness. A large \(\alpha\) gives regular, grid-like cells. A small \(\alpha\) lets cells follow colour edges closely.

Superpixels are usually a first step. To find whole objects, neighbouring superpixels with similar colour are merged.

The region adjacency graph (RAG) organises this:

• Vertices are regions (pixels or superpixels).
• Edges join regions that touch.
• Edge weights measure how different two regions are, for example the distance between their mean colours.

Iterative merging (slide 36) repeats one move:

1. Find the lowest-weight edge.
2. If its weight is below a threshold, merge its two regions into one.
3. Recompute the merged region's mean colour and the weights of its edges, and repeat.

Merging builds regions bottom-up. A cut works top-down: split the graph into two sets \(A\) and \(B\) by deleting edges.

Here edge weights measure similarity (large means alike). The cost of a split is the weight of the deleted edges:

$$\operatorname{cut}(A, B) = \sum_{i \in A,\ j \in B} w_{ij}.$$

The minimum cut has a flaw: cutting off one weakly connected node is cheap, so it favours tiny segments.

The normalised cut fixes this by dividing by each side's total connection to the whole graph, \(\operatorname{assoc}(A, V) = \sum_{i \in A} d_i\), where \(d_i\) is the weighted degree of node \(i\):⁶

$$\operatorname{Ncut}(A, B) = \frac{\operatorname{cut}(A, B)}{\operatorname{assoc}(A, V)} + \frac{\operatorname{cut}(A, B)}{\operatorname{assoc}(B, V)} .$$

A tiny \(A\) has a tiny \(\operatorname{assoc}(A, V)\), so its first term is large. Balanced splits along weak edges score best.

Clustering finds colour regions, and those regions often do not correspond to objects (slide 27).

The standard approach to instance segmentation (slide 42) has two stages:

1. Find a region of interest: a patch that is likely to contain one object.
2. Segment that patch into two classes: foreground (the object) and background.

Slide 43 gives two examples:

• Mask R-CNN uses a learned detector to propose boxes. For each box it predicts a class, a refined box, and a binary mask from a small fully-convolutional head. RoIAlign resamples each box to a fixed grid without rounding coordinates, so the mask lines up with the pixels.⁷ Section 4 of the notebook runs this model.
• GrabCut takes the box from a user. It models the colours of the foreground and of the background each as a Gaussian mixture, and labels pixels to maximise the likelihood plus a penalty for colour boundaries.⁸

GrabCut minimises an energy over the labels \(\alpha_p \in \{\text{fg}, \text{bg}\}\):

$$E(\boldsymbol\alpha) = \underbrace{\sum_{p} -\log p\big(\mathbf z_p \mid \alpha_p\big)}_{\text{colour likelihood}} \; + \; \lambda \underbrace{\sum_{(p,q)\ \text{neighbours}} [\alpha_p \ne \alpha_q]}_{\text{boundary length}} .$$

It alternates between refitting the two colour models and relabelling the pixels, and the energy falls at every round (slide 47).

Lecture 20 turns to methods that use labels. The first idea (slide 8) is to classify every pixel from an \(N \times N\) window around it.

Run naively, that is one CNN evaluation per pixel: \(480\,000\) evaluations for a small \(600 \times 800\) image (slide 9).

A fully-convolutional network (FCN) has no fully-connected layers. It produces a spatial map of class scores in one pass, so it classifies all windows in parallel and accepts any input size (slides 10–11).⁹

Slide 12 names two problems with a plain stack of convolutions:

• The receptive field grows only linearly. With \(3 \times 3\) kernels and stride 1, an output pixel sees \(2L + 1\) input pixels per side after \(L\) layers.
• Downsampling fixes that but loses resolution. Max pooling makes the receptive field grow geometrically, but the output map gets smaller.

The solution is an encoder–decoder (slide 13): the encoder downsamples, and the decoder upsamples back to \(H \times W\).

The decoder needs layers that make feature maps larger. Lecture 20 gives two.

Max unpooling

Each upsampling layer is paired with a max-pooling layer in the encoder (slide 15).¹⁰

1. When the encoder pools, it records where each maximum came from (the "switches").
2. When the decoder unpools, it writes each value back into that remembered position, and fills the rest of the block with zeros.

The simpler alternatives on slide 14 copy each value into the whole block (nearest neighbour) or always into the top-left corner (bed of nails).

Transposed convolution

A convolution with stride 2 is a learnable downsampler: the filter moves 2 pixels in the input for every 1 pixel in the output (slide 18).

A transposed convolution reverses the roles. The filter moves 2 pixels in the output for every 1 pixel in the input:

1. Multiply the whole kernel by one input value.
2. Stamp the result into the output at stride 2.
3. Where stamps overlap, add them (slide 20).
4. Trim the border to reach the target size (slide 21).

Some papers call this "deconvolution", "upconvolution" or "fractionally strided convolution".¹¹

Why the name "transposed"

Slides 23–24 write convolution as a matrix product. For kernel \(\vec x = [x, y, z]\), stride 2 and padding 1, convolving a length-4 signal \(\vec a\) (padded to length 6) is

$$\vec x * \vec a = X \vec a = \begin{bmatrix} x & y & z & 0 & 0 & 0 \\ 0 & 0 & x & y & z & 0 \end{bmatrix}\begin{bmatrix} 0 \\ a \\ b \\ c \\ d \\ 0 \end{bmatrix} = \begin{bmatrix} ay + bz \\ bx + cy + dz \end{bmatrix}.$$

The transposed convolution multiplies by the transpose of the same matrix. A length-2 input becomes a length-6 output, two of which are padding:

$$\vec x *^{T} \vec a = X^{T} \vec a = \begin{bmatrix} x & 0 \\ y & 0 \\ z & x \\ 0 & y \\ 0 & z \\ 0 & 0 \end{bmatrix}\begin{bmatrix} a \\ b \end{bmatrix} = \begin{bmatrix} ax \\ ay \\ az + bx \\ by \\ bz \\ 0 \end{bmatrix}.$$

Read down a column of \(X^{T}\): it is the kernel, shifted by the stride. That is exactly the "stamp and add" rule above.

The last layer. A \(1 \times 1\) convolution mixes channels at each pixel without looking at neighbours. With \(C\) output channels it turns every pixel's feature vector into \(C\) class scores (slide 26).

The loss. Training compares the predicted class distribution \(\hat{\mathbf y}_p\) with the true class at every pixel \(p\), and averages the cross-entropy over the \(N\) pixels (slide 27):

$$E = -\frac{1}{N}\sum_{p=1}^{N} \sum_{c=1}^{C} y_{p,c}\, \log \hat y_{p,c} .$$

U-Net was built for biomedical images.¹² It is an encoder–decoder with one extra idea: each decoder stage concatenates the upsampled map with the feature map from the matching encoder stage (slides 33–34).

The encoder copy carries fine spatial detail that pooling threw away, and the decoder supplies context. For cell segmentation, U-Net also weights the loss so that the thin gaps between touching cells count more.

The original network uses unpadded ("valid") \(3 \times 3\) convolutions, so each one shrinks the map by 2 pixels. The encoder maps are cropped before they are concatenated, which is the "copy and crop" arrow.

The Segment Anything Model (SAM) takes an image plus a prompt and returns a binary mask (slide 38).¹³ A prompt can be a point, a box, a mask or a text embedding.

Its encoder–decoder structure (slides 39–41) has three parts:

• Image encoder. A vision transformer pre-trained as a masked autoencoder¹⁴ turns the image into a \(64 \times 64\) grid of 256-dimensional vectors. It runs once per image.
• Prompt encoder. Points and boxes become positional embeddings, text goes through CLIP, and mask prompts go through a small convolution.
• Mask decoder. A few transformer layers let the prompt tokens and the image embedding attend to each other. The decoder then upsamples \(2\times\) with transposed convolutions, and outputs a mask logit at every pixel plus a predicted IoU score.

Ambiguity. A point on a person's shirt could mean the shirt or the whole person. SAM outputs three masks (for example subpart, part and whole), and is trained by back-propagating only the best of the three (slide 42).

Training. The mask loss combines focal loss and Dice loss, and the IoU head is trained with the mean squared error to the true IoU (slide 43). The SA-1B dataset has 11 million images and 1.1 billion masks, 99% generated automatically (slide 44).

Measuring overlap

Let \(P\) be the predicted pixels and \(G\) the ground-truth pixels. Two scores compare them:

$$\operatorname{IoU}(P, G) = \frac{|P \cap G|}{|P \cup G|}, \qquad \operatorname{Dice}(P, G) = \frac{2\,|P \cap G|}{|P| + |G|}.$$

Both are 1 for a perfect match and 0 for no overlap. Dice is always at least as large as IoU.

Focal loss down-weights pixels that are already easy. If \(p_t\) is the predicted probability of the true class at a pixel, it multiplies the cross-entropy by \((1 - p_t)^\gamma\):¹⁵

$$\operatorname{FL}(p_t) = -(1 - p_t)^{\gamma}\, \log p_t .$$

The Dice loss is \(1 - \operatorname{Dice}\), computed with soft probabilities.¹⁶ Because it only counts the object and the prediction, a huge background does not swamp it.

Lecture 20 ends with the trade-off for encoder–decoder methods (slide 47). They handle complex objects and backgrounds, and can take text as input or output. But their masks may miss precise boundaries, and they may miss background or hard-to-name objects.

The notebook starts by clustering the RGB colours of an aerial photo of Albert Park into \(K\) groups, then repainting each pixel with its cluster centre.

Cell numbers count every cell of worksheet11_solution-old-24.ipynb from 0, markdown included. Code is shown verbatim.

## For colab user, run
! pip install scikit-image==0.21.0
## After the update, You must restart the runtime in order to use newly installed versions.

import matplotlib.pyplot as plt
import numpy as np
import cv2
import tensorflow as tf
import logging

# Control the logging verbosity
# Suppress verbose Tensorflow logging
tf.get_logger().setLevel(logging.ERROR)
# Sets the AutoGraph verbosity level.
tf.autograph.set_verbosity(0)
# Verbosity is now 0

from sklearn.cluster import KMeans
from skimage.data import astronaut, coffee
from skimage.color import rgb2gray, label2rgb
from skimage.filters import sobel
from skimage import graph
from skimage.segmentation import felzenszwalb, slic, quickshift, watershed
from skimage.segmentation import mark_boundaries
from skimage.util import img_as_float
from sklearn.feature_extraction import image
from sklearn.cluster import spectral_clustering

# Read an image into BGR Format
!wget https://raw.githubusercontent.com/saraao/COMP90086_image/main/Albert_Park.jpg -q -O Albert_Park.jpg
img =  cv2.imread('Albert_Park.jpg')
# Convert to RGB
img = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)

print("Image shape:", img.shape)

plt.imshow(img) # as RGB Format
plt.axis("off")
plt.show()

Image shape: (800, 640, 3)

OpenCV reads images in BGR channel order, which is why cell 8 converts to RGB before plotting.

KMeans expects one row per sample, so cell 10 flattens the \(800 \times 640 \times 3\) array into \(512\,000\) rows of 3 colour values:

# Reshape the image to a 2D array of pixels and 3 color values
flat_image = img.reshape((-1,3))
# Convert to float
flat_image = np.float32(flat_image)
print(flat_image.shape)

(512000, 3)

K = 4
kmeans = KMeans(n_clusters=K)
kmeans.fit(flat_image)

kmeans.cluster_centers_

array([[125.390465, 162.35413 , 170.89093 ],
       [189.08794 , 198.64023 , 204.28804 ],
       [159.2623  , 154.95505 ,  88.992035],
       [ 60.432976,  83.3275  ,  75.4993  ]], dtype=float32)

# Convert to Integer format
centers = np.uint8(kmeans.cluster_centers_)
print(centers) # In RGB Format

plt.figure(0,figsize=(8,2))

i = 1
for each_col in centers:
    plt.subplot(1,4,i)
    plt.axis("off")
    i+=1

    # Color Swatch
    a = np.zeros((100,100,3),dtype='uint8')
    a[:,:,:] = each_col
    plt.imshow(a)

plt.show()

res = centers[kmeans.labels_]
result_image1 = res.reshape((img.shape))

plt.imshow(result_image1)
plt.title("Segmented Image when K = 4")
plt.axis("off")
plt.show()

# Use different K

# image 2
K = 10
kmeans = KMeans(n_clusters=K)
kmeans.fit(flat_image)
centers = np.uint8(kmeans.cluster_centers_)
res = centers[kmeans.labels_]
result_image2 = res.reshape((img.shape))

# image 3
K = 2
kmeans = KMeans(n_clusters=K)
kmeans.fit(flat_image)
centers = np.uint8(kmeans.cluster_centers_)
res = centers[kmeans.labels_]
result_image3 = res.reshape((img.shape))

# Plot them out
figure_size = 8
plt.figure(figsize=(figure_size,figure_size))
#original image
plt.subplot(2,2,1),plt.imshow(img)
plt.title('Original Image'), plt.xticks([]), plt.yticks([])
#image 2
plt.subplot(2,2,2),plt.imshow(result_image2)
plt.title('Segmented Image when K = 10'), plt.xticks([]), plt.yticks([])
#image 1
plt.subplot(2,2,3),plt.imshow(result_image1)
plt.title('Segmented Image when K = 4'), plt.xticks([]), plt.yticks([])
#image 3
plt.subplot(2,2,4),plt.imshow(result_image3)
plt.title('Segmented Image when K = 2'), plt.xticks([]), plt.yticks([])
plt.show()

Cell 18 is the key line. kmeans.labels_ holds one cluster index per pixel, and indexing centers with it replaces every pixel by its centre's colour. This is "colour quantisation".

Cell 23 sums it up: more detail needs more clusters. With \(K = 2\), the water and the grass merge. With \(K = 10\), the tree shadows and the road markings appear.

The clusters are colour groups, not objects. The road, the path and the brightest water share one cluster because their colours are close.

Cell 25 builds a five-dimensional feature vector for every pixel:

$$\mathbf f = \big(R,\ G,\ B,\ s\,x,\ s\,y\big), \qquad s = \texttt{xy\_scaling}.$$

scikit-learn's MeanShift uses one flat window of radius \(h\) for all five dimensions. Scaling the coordinates by \(s\) has the same effect as giving position its own bandwidth \(h / s\).

With \(h = 25\) and \(s = 0.005\), the spatial bandwidth is \(25 / 0.005 = 5000\) pixels, far larger than the \(256 \times 171\) image. So in the notebook's setting, position barely matters.

Running mean shift from all \(43\,776\) pixels would be slow. scikit-learn's implementation saves work in three steps:

1. Bin seeding. Round every feature vector to a grid of spacing \(h\), and start one window from each occupied grid cell. In the notebook's setting this leaves only a few hundred seeds.
2. Shift. Move each window to the mean of the points within distance \(h\) (a flat kernel), until it moves less than \(10^{-3} h\).
3. Remove duplicates. Sort the converged windows by how many points they contain. Keep one, and drop every other window within distance \(h\) of a kept one.

With cluster_all=True, every pixel then joins its nearest surviving centre.

from sklearn.cluster import MeanShift, estimate_bandwidth

# Read an image into BGR Format
!wget https://raw.githubusercontent.com/saraao/COMP90086_image/main/256kodim23.png -q -O 256kodim23.png
img =  cv2.imread('256kodim23.png')
# Convert to RGB
img = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)

# Show the image
plt.imshow(img)
plt.axis("off")
plt.show()

# Reshape the image to a 2D array of pixels and 3 color values
flat_image = img.reshape((-1,3))
# Convert to float
flat_image = np.float32(flat_image)

# Get x,y coordinate values of pixels
x, y = np.meshgrid(np.arange(img.shape[0]),np.arange(img.shape[1]))
x = x.reshape((-1,1))
y = y.reshape((-1,1))
xy = np.concatenate((x,y),axis=1)

# Since the sklearn function only allows one value for bandwidth, we will scale
# the (x,y) values relative to the colour values. This has the same effect as
# using different bandwidth values for space and colour.
xy_scaling = 0.005
xy = xy_scaling*xy

# Combine colour and x,y values
flat_image = np.concatenate((flat_image,xy),axis=1)

print(flat_image.shape)

(43776, 5)

# Estimate the bandwidth to use with the mean-shift algorithm.
bandwidth = estimate_bandwidth(flat_image,
                                 quantile=.04,
                                 n_jobs=-1, #use all processors your computer has
                                 n_samples=500) #The number of samples to use. Set a small value to reduce the time

bandwidth = 25

# Compute clustering with MeanShift
ms = MeanShift(bandwidth = bandwidth,
               n_jobs=-1,  #use all processors which your computer has
               bin_seeding=True, #True to speed up the algorithm because fewer seeds will be initialized
               cluster_all=True) #If true, then all points are clustered.
ms.fit(flat_image)

# Segmentation of the original image
centers = np.uint8(ms.cluster_centers_) # Convert centers to Integer format
res = centers[ms.labels_]
res = res[:,:3] # For visualisation, take just the color values (first 3 columns)
result_image = res.reshape((img.shape))

# Plot result
figure_size = 8
plt.figure(figsize=(figure_size,figure_size))
# original image
plt.subplot(1,2,1),plt.imshow(img)
plt.title('Original Image'), plt.xticks([]), plt.yticks([])
# final image
plt.subplot(1,2,2),plt.imshow(result_image)
plt.title(str(len(ms.cluster_centers_))+" colors,Mean-Shift"), plt.xticks([]), plt.yticks([])
plt.show()

# Compute clustering with MeanShift
# your code

# Segmentation of the original image
# your code

# Plot result
# your code

The figure below reruns cells 25 and 29 in your browser. It is a JavaScript port of scikit-learn's algorithm (bin seeding, flat kernel, duplicate removal, cluster_all=True).

On the same pixels, it finds the same number of clusters as scikit-learn 1.7 for every setting we tested: 14, 6, 67 and 32 clusters for \((h, s) = (25, 0.005), (40, 0.005), (15, 0.005), (25, 0.3)\).

Cell 31 compares four over-segmentation methods on the astronaut image, downsampled to \(256 \times 256\):

• Felzenszwalb–Huttenlocher. Treats pixels as a graph and merges greedily, using a scale parameter that favours larger regions.¹⁷
• Quickshift. An approximation of kernelised mean shift in the 5-D colour-and-position space.¹⁸ max_dist cuts its hierarchy at a chosen level.
• SLIC. Local k-means in \((l, a, b, x, y)\), as in A4.
• Compact watershed. Floods the Sobel gradient image from 250 markers. compactness penalises flooding far from a marker.¹⁹

# Reshape the astronaut image and convert the image to floating point format.
img = img_as_float(astronaut()[::2, ::2])

# Segments image using different methods
segments_fz = felzenszwalb(img, scale=100, sigma=0.5, min_size=50)
segments_slic = slic(img, n_segments=250, compactness=10, sigma=1,
                     start_label=1)
segments_quick = quickshift(img, kernel_size=3, max_dist=6, ratio=0.5)
gradient = sobel(rgb2gray(img))
segments_watershed = watershed(gradient, markers=250, compactness=0.001)

# Count the number of segments
print(f'Felzenszwalb number of segments: {len(np.unique(segments_fz))}')
print(f'SLIC number of segments: {len(np.unique(segments_slic))}')
print(f'Quickshift number of segments: {len(np.unique(segments_quick))}')
print(f'Watershed number of segments: {len(np.unique(segments_watershed))}')

# Plot them out
figure_size = 8
plt.figure(figsize=(figure_size,figure_size))
# Felzenszwalbs's method
# mark_boundaries allows us to see where we get the end of border of each segment
plt.subplot(2,2,1),plt.imshow(mark_boundaries(img, segments_fz))
plt.title("Felzenszwalbs's method"), plt.axis('off')
# SLIC
plt.subplot(2,2,2),plt.imshow(mark_boundaries(img, segments_slic))
plt.title('SLIC'), plt.axis('off')
# Quickshift
plt.subplot(2,2,3),plt.imshow(mark_boundaries(img, segments_quick))
plt.title('Quickshift'), plt.axis('off')
# Compact watershed
plt.subplot(2,2,4),plt.imshow(mark_boundaries(img, segments_watershed))
plt.title('Compact watershed'), plt.axis('off')
plt.show()

Felzenszwalb number of segments: 194
SLIC number of segments: 196
Quickshift number of segments: 695
Watershed number of segments: 256

Section 3 over-segments the coffee-cup image with SLIC (about 400 superpixels), builds a region adjacency graph, and then either merges regions (3.1) or cuts the graph (3.2).

def _weight_mean_color(g, src, dst, n):
    """Callback to handle merging nodes by recomputing mean color.

    The method expects that the mean color of `dst` is already computed.

    Parameters
    ----------
    g : RAG
        The graph under consideration.
    src, dst : int
        The vertices in `graph` to be merged.
    n : int
        A neighbor of `src` or `dst` or both.

    Returns
    -------
    data : dict
        A dictionary with the `"weight"` attribute set as the absolute
        difference of the mean color between node `dst` and `n`.
    """

    diff = g.nodes[dst]['mean color'] - g.nodes[n]['mean color']
    diff = np.linalg.norm(diff)
    return {'weight': diff}


def merge_mean_color(g, src, dst):
    """Callback called before merging two nodes of a mean color distance graph.

    This method computes the mean color of `dst`.

    Parameters
    ----------
    graph : RAG
        The graph under consideration.
    src, dst : int
        The vertices in `graph` to be merged.
    """
    g.nodes[dst]['total color'] += g.nodes[src]['total color']
    g.nodes[dst]['pixel count'] += g.nodes[src]['pixel count']
    g.nodes[dst]['mean color'] = (g.nodes[dst]['total color'] /
                                      g.nodes[dst]['pixel count'])

The two callbacks define the merge rule of A5:

• merge_mean_color runs just before two nodes merge. It adds up their total colour and pixel counts, so the merged region's mean is the pixel-weighted mean.
• _weight_mean_color sets the new weight between the merged node and each neighbour \(n\): the Euclidean distance between their mean RGB colours.

$$w(\text{dst}, n) = \big\lVert \bar{\mathbf c}_{\text{dst}} - \bar{\mathbf c}_{n} \big\rVert_2 .$$

# Load the Coffee cup image
img = coffee()

# Use the SLIC algorithm to segment the input image and obtain the region labels
labels = slic(img, compactness=30, n_segments=400, start_label=1)
# Compute the Region Adjacency Graph using mean colors.
g = graph.rag_mean_color(img, labels)
out = label2rgb(labels, img, kind='avg', bg_label=0)
out = np.uint8(out)

# Perform hierarchical merging of a RAG.
# Greedily merges the most similar pair of nodes until no edges lower than 'thresh' remain.
labels2 = graph.merge_hierarchical(labels, g,
                                   thresh=35, #Regions connected by an edge with weight smaller than thresh are merged.
                                   rag_copy=False,
                                   in_place_merge=True,
                                   merge_func=merge_mean_color,#This function is called before merging two nodes.
                                   weight_func=_weight_mean_color)#The function to compute the new weights of the nodes adjacent to the merged node.

# Return an RGB image where color-coded labels are painted over the image.
out2 = label2rgb(labels2, img, kind='avg', bg_label=0)
out2 = np.uint8(out2)

fig, ax = plt.subplots(nrows=2, sharex=True, sharey=True, figsize=(6, 8))

ax[0].imshow(out)
ax[1].imshow(out2)

for a in ax:
    a.axis('off')

plt.tight_layout()

merge_hierarchical keeps a priority queue of edges. It repeatedly pops the lowest-weight edge, merges its two regions if the weight is below thresh, and pushes the recomputed edges back. The live figure in Exercise 2 is a JavaScript port of this loop.

img = coffee()

labels = slic(img, compactness=30, n_segments=400,
                            start_label=1)
out = label2rgb(labels, img, kind='avg', bg_label=0)

g = graph.rag_mean_color(img, labels, mode='similarity')
labels_nc = graph.cut_normalized(labels, g)
out_nc = label2rgb(labels_nc, img, kind='avg', bg_label=0)

out = np.uint8(out)
out_nc = np.uint8(out_nc)

# Plot results
fig, ax = plt.subplots(nrows=2, sharex=True, sharey=True, figsize=(6, 8))

ax[0].imshow(out)
ax[1].imshow(out_nc)

for a in ax:
    a.axis('off')

plt.tight_layout()

Here mode='similarity' turns colour distances into similarities, \(w_{ij} = e^{-d_{ij}^2/\sigma}\) with \(\sigma = 255\), because a cut needs weights where large means alike (A6).

The threshold has the units of the edge weights: a distance between two mean RGB colours on the 0–255 scale.

The largest possible distance, from black to white, is

$$\lVert (255, 255, 255) - (0, 0, 0) \rVert_2 = 255\sqrt 3 \approx 441.7 .$$

So thresh=35 means: merge neighbours whose average colours differ by less than about 8% of the full range.

A segmentation can fail in two opposite ways:

• Over-segmentation. One object is split into many regions. Every boundary of the object is probably present, plus many false ones.
• Under-segmentation. Different objects share one region. Some true boundaries are lost, and later steps cannot recover them.

A merging threshold moves along this spectrum: \(\tau = 0\) keeps every superpixel, and a very large \(\tau\) merges the whole image into one region.

Section 4 runs pre-trained models from Detectron2,²⁰ Meta's detection library, on a COCO validation image.²¹ It needs a GPU runtime on Colab.

#@title Check gpu
!nvidia-smi

#@title Install detectron2
!python -m pip install 'git+https://github.com/facebookresearch/detectron2.git'

#@title Check PyTorch, CUDA and Detectron2 Versions
import torch, detectron2
!nvcc --version
TORCH_VERSION = ".".join(torch.__version__.split(".")[:2])
CUDA_VERSION = torch.__version__.split("+")[-1]
print("torch: ", TORCH_VERSION, "; cuda: ", CUDA_VERSION)
print("detectron2:", detectron2.__version__)

# You may need to restart your runtime prior to this, to let your installation take effect
# Some basic setup:
# Setup detectron2 logger
import detectron2
from detectron2.utils.logger import setup_logger
setup_logger()

# import some common libraries
import numpy as np
import os, json, cv2, random
from google.colab.patches import cv2_imshow

# import some common detectron2 utilities
from detectron2 import model_zoo
from detectron2.engine import DefaultPredictor
from detectron2.config import get_cfg
from detectron2.utils.visualizer import Visualizer
from detectron2.data import MetadataCatalog, DatasetCatalog

!wget http://images.cocodataset.org/val2017/000000439715.jpg -q -O input.jpg
im = cv2.imread("./input.jpg")
cv2_imshow(im)

Cell 49 loads a Mask R-CNN with a ResNet-50 feature-pyramid backbone, trained on COCO's 80 categories (A7). Detections with a score below 0.5 are discarded.

cfg = get_cfg()
# add project-specific config (e.g., TensorMask) here if you're not running a model in detectron2's core library
cfg.merge_from_file(model_zoo.get_config_file("COCO-InstanceSegmentation/mask_rcnn_R_50_FPN_3x.yaml"))
cfg.MODEL.ROI_HEADS.SCORE_THRESH_TEST = 0.5  # set threshold for this model
# Find a model from detectron2's model zoo. You can use the https://dl.fbaipublicfiles... url as well
cfg.MODEL.WEIGHTS = model_zoo.get_checkpoint_url("COCO-InstanceSegmentation/mask_rcnn_R_50_FPN_3x.yaml")
predictor = DefaultPredictor(cfg)
outputs = predictor(im)

# look at the outputs. See https://detectron2.readthedocs.io/tutorials/models.html#model-output-format for specification
print(outputs["instances"].pred_classes)
print(outputs["instances"].pred_boxes)

tensor([17,  0,  0,  0,  0,  0,  0,  0, 25,  0, 25, 25,  0,  0, 24],
       device='cuda:0')
Boxes(tensor([[126.6035, 244.8977, 459.8291, 480.0000],
        [251.1083, 157.8127, 338.9731, 413.6379],
        [114.8496, 268.6864, 148.2352, 398.8111],
        [  0.8217, 281.0327,  78.6072, 478.4210],
        [ 49.3953, 274.1229,  80.1545, 342.9808],
        [561.2248, 271.5816, 596.2755, 385.2552],
        [385.9072, 270.3125, 413.7130, 304.0397],
        [515.9295, 278.3743, 562.2792, 389.3803],
        [335.2409, 251.9167, 414.7491, 275.9375],
        [350.9300, 269.2060, 386.0984, 297.9081],
        [331.6292, 230.9996, 393.2759, 257.2009],
        [510.7349, 263.2656, 570.9865, 295.9194],
        [409.0841, 271.8646, 460.5582, 356.8722],
        [506.8767, 283.3257, 529.9403, 324.0392],
        [594.5663, 283.4820, 609.0577, 311.4124]], device='cuda:0'))

Boxes are \((x_1, y_1, x_2, y_2)\) in pixels. Classes are zero-based COCO indices without the background class: 0 is person, 17 horse, 24 backpack and 25 umbrella.

# We can use `Visualizer` to draw the predictions on the image.
v = Visualizer(im[:, :, ::-1], MetadataCatalog.get(cfg.DATASETS.TRAIN[0]), scale=1.2)
out = v.draw_instance_predictions(outputs["instances"].to("cpu"))
cv2_imshow(out.get_image()[:, :, ::-1])

A panoptic model labels every pixel. Things (person, horse, umbrella) get instance masks, and stuff (tree, fence, dirt) gets one region per class (A1).

# Inference with a panoptic segmentation model
cfg = get_cfg()
cfg.merge_from_file(model_zoo.get_config_file("COCO-PanopticSegmentation/panoptic_fpn_R_101_3x.yaml"))
cfg.MODEL.WEIGHTS = model_zoo.get_checkpoint_url("COCO-PanopticSegmentation/panoptic_fpn_R_101_3x.yaml")
predictor = DefaultPredictor(cfg)
panoptic_seg, segments_info = predictor(im)["panoptic_seg"]
v = Visualizer(im[:, :, ::-1], MetadataCatalog.get(cfg.DATASETS.TRAIN[0]), scale=1.2)
out = v.draw_panoptic_seg_predictions(panoptic_seg.to("cpu"), segments_info)
cv2_imshow(out.get_image()[:, :, ::-1])

# Inference with a panoptic segmentation model
cfg = get_cfg()
cfg.merge_from_file(model_zoo.get_config_file("COCO-PanopticSegmentation/panoptic_fpn_R_101_3x.yaml"))
cfg.MODEL.WEIGHTS = model_zoo.get_checkpoint_url("COCO-PanopticSegmentation/panoptic_fpn_R_101_3x.yaml")
predictor = DefaultPredictor(cfg)
panoptic_seg, segments_info = predictor(im)["panoptic_seg"]
v = Visualizer(im[:, :, ::-1], MetadataCatalog.get(cfg.DATASETS.TRAIN[0]), scale=1.2)
out = v.draw_panoptic_seg_predictions(panoptic_seg.to("cpu"), segments_info)
cv2_imshow(out.get_image()[:, :, ::-1])

Mask R-CNN segments only what its detector proposes and its classifier accepts:

1. The region proposal network suggests boxes that look like COCO objects.
2. Each box is classified into one of 80 classes or background, with a score.
3. Only boxes with a score above SCORE_THRESH_TEST = 0.5 survive, and only those get a mask.

Anything outside the 80 classes can only be "background", or be forced into the nearest class. A low score means no mask at all, even when the object is plainly visible.

Detectors produce many overlapping boxes for one object. Non-maximum suppression keeps the highest-scoring box and deletes any other box of the same class whose IoU with it exceeds a threshold (0.5 in Detectron2's test configuration).

For two axis-aligned boxes, the intersection is itself a box:

$$|P \cap G| = \max\!\big(0,\ \min(x_2, x_2') - \max(x_1, x_1')\big)\cdot \max\!\big(0,\ \min(y_2, y_2') - \max(y_1, y_1')\big).$$

In a crowd, two different people can overlap by more than the threshold, and one of them is suppressed.

Mean shift on colour (with no position features) groups pixels by colour alone. Two pixels anywhere in the image join the same cluster if their colours climb to the same peak.

So a cluster is a set of colours, not a connected object. An object made of two very different colours, like a black-and-white penguin, is split between two clusters, and each of those clusters also contains every background pixel of the same colour.

Superpixels are spatially compact, and merging only joins adjacent regions. So separate black patches never merge with each other unless they touch, which is an improvement.

But the merge test is still colour similarity. A white belly touching white sand looks like one region.

Points worth raising in the tutorial, found while porting the notebook:

1. Cell 25 builds the coordinates in the wrong order. np.meshgrid(np.arange(H), np.arange(W)) has shape \((W, H)\), so after flattening the \((x, y)\) values do not belong to the pixels they are attached to (the image is \(171 \times 256\), not square). Use y, x = np.mgrid[0:H, 0:W]. It is invisible at xy_scaling = 0.005; at 1 the clusters turn into horizontal stripes (Exercise 1 figure, with false colours).
2. With xy_scaling = 0.005 the spatial term does almost nothing: the coordinates span at most about 1.3 units against a bandwidth of 25.
3. Cell 28's bandwidth estimate is overwritten by bandwidth = 25 in cell 29.
4. Cells 52 and 54 are identical. The markdown in cell 53 promises "other types of models", but the panoptic model simply runs again.
5. Cell 12 creates KMeans without n_init or random_state, so it prints a FutureWarning and the cluster order (and the swatch order in cell 16) changes from run to run. Cell 16 also hard-codes 4 subplots.
6. The docstring of _weight_mean_color says "absolute difference", but the code takes the Euclidean norm of the RGB difference.
7. Small typos in cell 30 and 33: "Felsenszwalb" for Felzenszwalb, and the label2rgb link points to lab2rgb.
8. cut_normalized is not deterministic (random eigensolver start), so cell 37 can differ between runs.
9. Exercises 2 and 3 have no solution cells (cells 39 and 57 are empty). The answers on this page are ours.