Spatial Clustering of Archaeological Sites

The Problem

• Archaeological surveys record the locations of sites (settlements, hearths, burial mounds) as 2-D coordinates on a map.
• Sites are rarely distributed uniformly: concentrations may represent settlements around a water source, nodes in a trade network, or ritual landscapes.
• Spatial clustering groups sites by proximity without requiring the analyst to specify the number of groups in advance.
• The approach here:
• Generate a synthetic set of site locations in which clusters at known centers are surrounded by random background noise.
• Divide the survey region into a regular grid and count sites per cell.
• Mark cells whose count meets a density threshold as "hot."
• Find connected components of hot cells using depth-first search (DFS) with 8-connectivity.
• Assign each site the label of the component containing its cell.
• Compute the centroid of each cluster.

Dividing the Region into a Grid

• The survey region spans $[\text{region_min}, \text{region_max}]$ in both easting and northing.
• Divide each axis into grid_size equal-width cells of width $w = (\text{region_max} - \text{region_min}) / \text{grid_size}$.
• The cell index along one axis for a coordinate $v$ is:

$$\text{cell} = \left\lfloor \frac{v - \text{region_min}}{w} \right\rfloor$$

• Sites exactly on the upper boundary would map to index grid_size; clip to grid_size - 1 to keep them inside the array.
• Count the number of sites in each cell; the result is a 2D integer array of shape (grid_size, grid_size).

Identifying Concentration Cells

• A cell with count >= HOT_THRESHOLD is called a "hot" cell; it represents a local concentration of sites.
• Cells below the threshold are background: isolated or sparse finds that do not indicate a settlement pattern.
• The Boolean mask of hot cells is straightforward:

• Raising HOT_THRESHOLD requires denser concentrations before a cell qualifies, which can split a real cluster into fragments or eliminate it entirely.
• Lowering HOT_THRESHOLD to 1 makes every occupied cell hot, merging isolated noise points into spurious clusters.

Finding Clusters with Depth-First Search

• Two hot cells belong to the same cluster if they are connected through a chain of adjacent hot cells.
• 8-connectivity counts all eight neighbors (four cardinal, four diagonal), so clusters that touch only at corners are still merged.
• The same iterative DFS used in the porous lesson labels connected components; here it assigns cluster IDs instead of testing for top-to-bottom percolation.

• labels is initialized to -1; only hot cells receive a non-negative label.
• Each cell is pushed onto the stack at most once (when first labelled), so the algorithm runs in $O(\text{grid_size}^2)$ time.

Assigning Sites to Clusters

• Each site is mapped to a grid cell using the same floor-division formula as build_grid.
• The site receives the cluster label of its cell; sites in background cells receive -1 (noise).

• Sites in hot cells that belong to a small disconnected cluster receive a distinct positive label, not noise.
• Two sites in adjacent cells may receive different labels if one cell is hot and the other is not, even if the sites are physically close.

Cluster Centroids

• The centroid of a cluster is the mean easting and mean northing of all sites assigned to that cluster.
• Centroids summarize each cluster's geographic location and can be compared with known landscape features.

• Noise sites (label -1) are excluded from all centroid calculations.
• The centroid of a cluster recovered from Gaussian-distributed data should be close to the planted center; the offset is the estimation error.

Testing

Grid shape and total count

• build_grid must return an array of shape (grid_size, grid_size).
• The sum of all cell counts must equal the total number of sites, because every site lands in exactly one cell.

Hot-cell mask

• find_hot_cells must mark cells with count >= threshold as True and all other cells as False.

Two disconnected blocks become two clusters

• A 3x3 hot mask with True at (0, 0) and True at (2, 2) must receive two distinct non-negative cluster IDs, because those cells are not 8-connected.

At least two clusters on synthetic data

• Running the full pipeline on the synthetic data must produce at least two clusters; the four planted concentrations are far enough apart to remain separate under the default parameters.

One centroid per label

• cluster_centroids must return exactly one entry for every unique non-negative label in the input, and the centroid coordinates must equal the mean of the corresponding sites.

import numpy as np

from generate_archcluster import (
    make_sites,
    REGION_MIN,
    REGION_MAX,
)
from archcluster import (
    build_grid,
    find_hot_cells,
    find_clusters,
    assign_labels,
    cluster_centroids,
    GRID_SIZE,
    HOT_THRESHOLD,
)


def test_build_grid_shape():
    # Grid must have the requested number of rows and columns.
    coords = np.array([[10.0, 20.0], [50.0, 50.0], [90.0, 80.0]])
    grid = build_grid(coords, REGION_MIN, REGION_MAX, GRID_SIZE)
    assert grid.shape == (GRID_SIZE, GRID_SIZE)


def test_build_grid_sum():
    # Every site must land in exactly one cell, so the grid sum equals site count.
    df = make_sites()
    coords = df.select(["easting", "northing"]).to_numpy()
    grid = build_grid(coords, REGION_MIN, REGION_MAX, GRID_SIZE)
    assert grid.sum() == len(coords)


def test_find_hot_cells_threshold():
    # Cells with count >= threshold must be True; cells below must be False.
    grid = np.array([[0, 1, 2], [3, 4, 0], [1, 2, 5]])
    hot = find_hot_cells(grid, threshold=2)
    expected = grid >= 2
    assert np.array_equal(hot, expected)


def test_find_clusters_two_disconnected_blocks():
    # Two separated 1-cell hot regions must receive distinct non-negative IDs.
    hot = np.zeros((3, 3), dtype=bool)
    hot[0, 0] = True   # top-left block
    hot[2, 2] = True   # bottom-right block (not 8-connected to top-left)
    labels = find_clusters(hot)
    lbl_a = labels[0, 0]
    lbl_b = labels[2, 2]
    assert lbl_a >= 0
    assert lbl_b >= 0
    assert lbl_a != lbl_b
    # Cold cells stay -1.
    assert labels[0, 1] == -1
    assert labels[1, 0] == -1


def test_find_clusters_8connectivity():
    # Diagonal neighbors must be merged into one cluster.
    hot = np.zeros((3, 3), dtype=bool)
    hot[0, 0] = True
    hot[1, 1] = True   # diagonal from (0,0)
    labels = find_clusters(hot)
    assert labels[0, 0] == labels[1, 1]
    assert labels[0, 0] >= 0


def test_assign_labels_produces_two_or_more_clusters():
    # On the standard synthetic data the pipeline must find at least two clusters.
    df = make_sites()
    coords = df.select(["easting", "northing"]).to_numpy()
    grid = build_grid(coords, REGION_MIN, REGION_MAX, GRID_SIZE)
    hot = find_hot_cells(grid, HOT_THRESHOLD)
    cell_labels = find_clusters(hot)
    labels = assign_labels(coords, cell_labels, REGION_MIN, REGION_MAX, GRID_SIZE)
    n_clusters = len(set(labels[labels >= 0]))
    assert n_clusters >= 2


def test_cluster_centroids_one_per_label():
    # cluster_centroids must return exactly one entry for every unique non-negative label.
    coords = np.array([
        [10.0, 10.0],
        [11.0, 10.0],
        [50.0, 50.0],
        [51.0, 50.0],
        [90.0, 90.0],  # noise
    ])
    labels = np.array([0, 0, 1, 1, -1])
    centroids = cluster_centroids(coords, labels)
    assert set(centroids.keys()) == {0, 1}
    cx0, cy0 = centroids[0]
    assert abs(cx0 - 10.5) < 1e-9
    assert abs(cy0 - 10.0) < 1e-9

DBSCAN handles clusters of irregular shape and density without discretizing space; it requires a more advanced algorithms course.