Run a Nearest Neighbour Analysis in PyQGIS
Looking at a map of points, almost everyone sees clusters. Human eyes find patterns in randomness, which is why the first question in point pattern analysis is a statistical one: are these points closer together than they would be by chance? The average nearest neighbour index answers it with a single number. It compares the mean distance from each point to its nearest neighbour with the mean distance expected if the same number of points were scattered at random over the same area.
This recipe belongs to Spatial Statistics & Pattern Analysis. It runs the nearest neighbour analysis algorithm, interprets its outputs, shows why the study area matters more than any other choice, and computes the per-point nearest distances that the summary is built from.
Prerequisites
- QGIS 3.34 LTR or newer, or the QGIS 4 series.
- A point layer in a projected CRS; distances in degrees are meaningless here.
- A clear definition of the study area — the region within which points could have occurred. This is not optional, as the section on study area explains.
Run the algorithm
The Processing algorithm takes a point layer and returns the observed and expected mean distances, the index, the number of points and a z-score, plus an HTML report.
import processing
from qgis.core import QgsProject
cafes = QgsProject.instance().mapLayersByName("cafes")[0]
result = processing.run("qgis:nearestneighbouranalysis", {
"INPUT": cafes,
"OUTPUT_HTML_FILE": "/data/results/cafes_nn.html",
})
print(f"points: {result['POINT_COUNT']}")
print(f"observed: {result['OBSERVED_MD']:.1f} m")
print(f"expected: {result['EXPECTED_MD']:.1f} m")
print(f"index: {result['NN_INDEX']:.3f}")
print(f"z-score: {result['Z_SCORE']:.2f}")
Breakdown: OBSERVED_MD is the mean of each point's distance to its nearest neighbour. EXPECTED_MD is the mean that complete spatial randomness would give for this many points in this area, computed as 0.5 divided by the square root of the density. NN_INDEX is their ratio: below 1 means points are closer than random (clustering), above 1 means further apart (dispersion). The z-score says how surprising the observed value is under randomness: beyond about ±1.96 the pattern differs from random at the 5 % level. The HTML file holds the same numbers in a form you can attach to a report.
Read the result correctly
The index and z-score answer different questions, and both are needed.
The index is an effect size: 0.4 is strong clustering, 0.9 is mild. The z-score is evidence: with thousands of points even a tiny departure from randomness gives a large z-score, while with twenty points a substantial-looking index can be indistinguishable from chance. Report both, and do not over-interpret a significant z-score attached to an index close to 1 — it means "detectably not random", not "strongly clustered".
Choose the study area deliberately
The expected distance depends on density, and density depends on the area you divide by. The algorithm, by default, uses the area of the points' bounding box. That is rarely the right area, and the choice can flip the conclusion.
import math
from qgis.core import QgsFeatureRequest
def nn_with_area(points, area_m2):
"""Clark–Evans index for a given study area."""
from qgis.core import QgsSpatialIndex
geoms = {f.id(): f.geometry() for f in points.getFeatures()}
index = QgsSpatialIndex(points.getFeatures())
dists = []
for fid, g in geoms.items():
nearest = [n for n in index.nearestNeighbor(g.asPoint(), 2) if n != fid][0]
dists.append(g.distance(geoms[nearest]))
n = len(dists)
observed = sum(dists) / n
expected = 0.5 / math.sqrt(n / area_m2)
se = 0.26136 / math.sqrt(n * n / area_m2)
return observed / expected, (observed - expected) / se
city = QgsProject.instance().mapLayersByName("city_boundary")[0]
city_area = sum(f.geometry().area() for f in city.getFeatures())
bbox_area = cafes.extent().area()
print("bounding box:", nn_with_area(cafes, bbox_area))
print("city boundary:", nn_with_area(cafes, city_area))
Breakdown: This reimplements the Clark–Evans statistic so the area can be chosen. With the bounding box, the expected distance assumes points could have appeared anywhere inside a rectangle around them; with the city boundary, anywhere in the city. If cafés are concentrated in the centre of a large city, the city area produces a much larger expected distance and a stronger clustering result than the bounding box. Neither is wrong — they answer different questions — but the area must match the question you are asking. Water bodies, parks and other places where points cannot occur should be removed from the area for an honest test.
Inspect distances point by point
The summary hides where clustering happens. Computing each point's nearest-neighbour distance as a field shows it on the map and allows a histogram.
import processing
nearest = processing.run("native:joinbynearest", {
"INPUT": cafes, "INPUT_2": cafes, "FIELDS_TO_COPY": ["name"],
"DISCARD_NONMATCHING": False, "PREFIX": "nn_", "NEIGHBORS": 1,
"MAX_DISTANCE": None, "OUTPUT": "memory:cafes_nn"})["OUTPUT"]
dists = sorted(f["distance"] for f in nearest.getFeatures())
median = dists[len(dists) // 2]
print(f"median {median:.0f} m, 90th percentile {dists[int(len(dists) * 0.9)]:.0f} m")
QgsProject.instance().addMapLayer(nearest)
Breakdown: Joining a layer to itself by nearest feature excludes each feature's match with itself, so distance holds the true nearest-neighbour distance and nn_name names the neighbour. Styling by that field with a graduated renderer shows dense cores and isolated outliers at a glance. The median is often a better summary than the mean, because a handful of isolated points can pull the mean far to the right. Joining attributes by nearest neighbour covers the algorithm's other options.
Know the limits of the test
Average nearest neighbour is a quick first test, not a complete analysis. It looks at one scale only — the distance to the very nearest point — so a pattern that is clustered at 200 m and regular at 2 km reports only the first. It is sensitive to edge effects: points near the study area boundary have their true nearest neighbour outside it, inflating observed distances. And it says nothing about where the clusters are.
For multiple scales, compute nearest distances to the kth neighbour or use Ripley's K function from a Python library. For where clusters are, use DBSCAN or k-means clustering or a hot spot analysis. For a visual density surface, kernel density estimation.
Correct for edge effects with a guard zone
Points near the edge of the study area are disadvantaged: their true nearest neighbour may lie just outside the boundary, in data you did not include, so their measured nearest distance is too long. With many points near the edge, this biases the index towards dispersion. A guard zone fixes it: measure distances for points inside an inner area only, but let them find neighbours anywhere in the full dataset.
GUARD = 300 # metres, about the typical nearest distance or larger
inner = next(city.getFeatures()).geometry().buffer(-GUARD, 8)
all_geoms = {f.id(): f.geometry() for f in cafes.getFeatures()}
index = QgsSpatialIndex(cafes.getFeatures())
dists = []
for fid, g in all_geoms.items():
if not inner.contains(g):
continue # edge point: may be a neighbour, not a subject
nearest = [n for n in index.nearestNeighbor(g.asPoint(), 2) if n != fid][0]
dists.append(g.distance(all_geoms[nearest]))
n_all = len(all_geoms)
observed = sum(dists) / len(dists)
expected = 0.5 / math.sqrt(n_all / city_area)
print(f"guarded: {len(dists)} subjects, index {observed / expected:.3f}")
Breakdown: Shrinking the boundary inward by the guard distance with a negative buffer defines the subjects; every point, including those in the guard zone, remains a candidate neighbour through the spatial index. The expected distance still uses the density of the whole area, which is what the subjects experience. Choose a guard width at least as large as typical nearest distances — the median from the per-point section is a good starting value. If the guarded index differs noticeably from the unguarded one, edge effects were distorting the result and the guarded value is the one to report.
Compare groups and time periods
The index becomes more useful when it is compared: cafés against restaurants in the same city, this year's burglaries against last year's, one district against another.
for category in ["cafe", "restaurant", "bar"]:
subset = cafes.materialize(QgsFeatureRequest().setFilterExpression(
f"\"amenity\" = '{category}'"))
idx, z = nn_with_area(subset, city_area)
print(f"{category:<12} n={subset.featureCount():>5} index={idx:.2f} z={z:.1f}")
Breakdown: Using the same study area for every group makes the indexes comparable; changing the area between groups invalidates the comparison. Group sizes matter too: a small group has a noisier index, so the z-score is the fairer guide to whether a difference is real. Materialising each subset is simplest for a handful of categories; for many, filter inside nn_with_area instead.
QGIS version compatibility
qgis:nearestneighbouranalysis and native:joinbynearest are present on QGIS 3.34 LTR, 3.40 LTR and QGIS 4 with the same outputs. QgsSpatialIndex.nearestNeighbor with a point argument works on all; in QGIS 4 it also accepts a QgsGeometry.
Troubleshooting
- The observed distance is tiny or zero. The layer contains duplicate points; remove them first with duplicate detection.
- Distances look like fractions. The layer is in degrees; reproject to a metric CRS.
- Results swing wildly between runs on subsets. Very few points; the test needs dozens at least.
- Every pattern comes out clustered. The study area is much larger than the region points can occur in; use a tighter, honest area.
Conclusion
Run qgis:nearestneighbouranalysis for a first answer, read the index as effect size and the z-score as evidence, recompute with a study area that matches your question, inspect per-point distances on the map, and move to clustering or hot spot methods when you need to know where the pattern is.
Frequently Asked Questions
What is a good minimum number of points? At least thirty for a stable result; fewer and the z-score is unreliable.
Does the test work for polygons? Use their centroids, but be aware that large polygons make centroid distances misleading.
Is a dispersed result interesting? Yes — it suggests competition or regulation, such as shops avoiding each other or planning rules spacing facilities.
Can I weight points? Not in this test. Use kernel density or hot spot analysis with a weight field.