Spatial Weights Warn About Islands or Moran's I Returns nan

Problem statement

Building spatial weights produces a warning you were not expecting:

UserWarning: The weights matrix is not fully connected:
 There are 3 disconnected components.
 There are 2 islands with ids: 64, 65.

Or the statistic itself comes back as nothing at all:

moran = Moran(zones["rate"].to_numpy(), w, permutations=999)
print(f"I = {moran.I}  p = {moran.p_sim}")
I = nan  p = 0.001

Look carefully at that second line. I is nan and the p-value is 0.001 β€” a normal, significant-looking number. Any pipeline that checks the p-value and reports "significant clustering" will sail straight past this.

These two problems look related and are not. Islands produce a warning and a quietly reduced analysis. nan comes from somewhere else entirely.

Quick answer

import numpy as np

values = zones["rate"].to_numpy()
print(f"islands       {len(w.islands)}  {w.islands[:5]}")
print(f"components    {w.n_components}")
print(f"NaN values    {np.isnan(values).sum()}")
print(f"inf values    {np.isinf(values).sum()}")
print(f"variance      {np.nanvar(values):.6f}")
islands       2  [64, 65]
components    3
NaN values    1
inf values    0
variance      12.481744
Symptom Cause Fix
"not fully connected" warning disconnected units drop them, or use KNN weights
islands with ids: … units with no neighbours same
I = nan a single NaN or inf in the values drop or impute before computing
I = nan zero variance (all values identical) there is nothing to correlate
p_sim normal while I is nan permutation p is computed regardless check np.isfinite(moran.I) explicitly

Islands do not cause nan. A dataset with two islands returned I = 0.7877 perfectly happily β€” the islands were simply given zero weight and excluded from the calculation without appearing anywhere in the result.

Two independent problems: islands producing a connectivity warning and silent exclusion, and a single NaN value producing I equals nan with a normal p-value.
Two unrelated failures that get confused because they surface at the same line of code.

Step-by-step solution

1. Find the nan first β€” it is the one that voids the result

values = zones["rate"].to_numpy()
bad = ~np.isfinite(values)
if bad.any():
    print(f"{bad.sum()} non-finite values at index {zones.index[bad].tolist()[:5]}")
    print(zones.loc[bad, ["name", "count", "population", "rate"]].head())
1 non-finite values at index [10]
      name  count  population  rate
10  Zone 10      3           0   inf

A single division by zero. Moran's I is computed from sums over all units, so one non-finite value poisons the entire statistic β€” nan propagates through every term.

The commonest sources, in order:

  • a rate with a zero denominator β€” a zone with no population
  • a failed join leaving NaN in the value column
  • suppressed values in official statistics, often -99999 rather than NaN, which do not raise but do wreck the mean

2. Check that zero variance is not the cause

if np.nanvar(values) == 0:
    print("all values identical β€” Moran's I is undefined")

Moran's I divides by the variance of the values. If every unit has the same value there is no variation to correlate spatially, and the result is 0/0. This happens more often than expected β€” a filtered subset where every row has the same category, or a count column that is zero everywhere.

3. Then deal with the islands

print(f"{len(w.islands)} islands of {w.n} units")
print(f"{w.n_components} disconnected components")
print(zones.loc[w.islands, ["name"]].head())
2 islands of 66 units
3 disconnected components
       name
64  Island A
65  Island B

An island is a unit that shares no boundary with any other β€” a literal island, an exclave, or an artefact of geometry that does not quite touch. Contiguity weights give it no neighbours, so it contributes nothing to the statistic and is silently dropped.

That is not always wrong. But it means your reported n is larger than the number of units the statistic actually used, and nothing in the output says so.

4. Choose a deliberate response to the islands

# a) drop them, and say so
connected = zones.drop(index=w.islands)
w_dropped = Queen.from_dataframe(connected, use_index=True)
w_dropped.transform = "r"
print(f"dropped {len(w.islands)}: {w_dropped.n} units, {len(w_dropped.islands)} islands")

# b) use k-nearest neighbours, which guarantees every unit has neighbours
w_knn = KNN.from_dataframe(zones, k=4)
w_knn.transform = "r"
print(f"KNN k=4: {w_knn.n} units, {len(w_knn.islands)} islands")
dropped 2: 64 units, 0 islands
KNN k=4: 66 units, 0 islands

Both are defensible, and they answer slightly different questions:

  Queen with islands   I = 0.7877   (islands silently excluded)
  Queen, islands dropped  I = 0.7934   (n = 64, stated)
  KNN k=4              I = 0.7846   (n = 66, islands connected by distance)

The three agree closely here, which is the reassuring case. When they disagree, the islands were doing real work and the choice needs justifying.

5. Assert the result is finite

moran = Moran(values, w, permutations=999)
if not np.isfinite(moran.I):
    raise ValueError("Moran's I is not finite β€” check for NaN, inf or zero variance")
print(f"I = {moran.I:+.4f}  p = {moran.p_sim:.4f}  n = {w.n - len(w.islands)}")

Reporting n as the number of units with neighbours, not the number of rows, is what makes the result honest.

A grid of connected cells plus two detached island cells that receive zero weight and drop out of the statistic without appearing in the output.
Islands do not error. They get zero weight, contribute nothing, and the reported n still counts them.

Code examples

Example 1 β€” a weights builder that refuses to hide anything

import geopandas as gpd
import numpy as np
import warnings
from esda.moran import Moran
from libpysal.weights import KNN, Queen, Rook


def build_weights(zones, scheme="queen", *, k=4, transform="r", islands="report"):
    """Build spatial weights and report exactly what happened to disconnected units."""
    builders = {
        "queen": lambda g: Queen.from_dataframe(g, use_index=True),
        "rook": lambda g: Rook.from_dataframe(g, use_index=True),
        "knn": lambda g: KNN.from_dataframe(g, k=k),
    }

    with warnings.catch_warnings():
        warnings.simplefilter("ignore")          # we report this ourselves, in detail
        w = builders[scheme](zones)

    if w.islands and islands == "drop":
        kept = zones.drop(index=w.islands)
        print(f"dropped {len(w.islands)} island(s): {list(w.islands)[:5]}")
        with warnings.catch_warnings():
            warnings.simplefilter("ignore")
            w = builders[scheme](kept)
        zones = kept
    elif w.islands:
        print(f"{len(w.islands)} island(s) will get zero weight and be excluded "
              f"from the statistic: {list(w.islands)[:5]}")

    w.transform = transform
    print(f"{scheme}: n={w.n}, {w.n_components} component(s), "
          f"mean neighbours {w.mean_neighbors:.2f}, "
          f"min {min(w.cardinalities.values())}, max {max(w.cardinalities.values())}")
    return w, zones


def safe_moran(zones, column, w, *, permutations=999):
    values = zones[column].to_numpy(dtype="float64")

    bad = ~np.isfinite(values)
    if bad.any():
        raise ValueError(
            f"{bad.sum()} non-finite value(s) in {column!r} at index "
            f"{zones.index[bad].tolist()[:5]} β€” Moran's I would return nan"
        )
    if np.var(values) == 0:
        raise ValueError(f"{column!r} has zero variance β€” Moran's I is undefined")

    moran = Moran(values, w, permutations=permutations)
    if not np.isfinite(moran.I):
        raise ValueError(f"Moran's I is {moran.I} despite finite inputs")

    effective = w.n - len(w.islands)
    print(f"I = {moran.I:+.4f}  E[I] = {moran.EI:+.4f}  z = {moran.z_sim:+.2f}  "
          f"p = {moran.p_sim:.4f}  (n = {effective})")
    return moran


w, zones_used = build_weights(zones, "queen")
moran = safe_moran(zones_used, "rate", w)
2 island(s) will get zero weight and be excluded from the statistic: [64, 65]
queen: n=66, 3 component(s), mean neighbours 6.36, min 0, max 8
ValueError: 1 non-finite value(s) in 'rate' at index [10] β€” Moran's I would return nan

The exception names the column, the count and the offending index. Compare that with I = nan, p = 0.001 slipping through a pipeline unnoticed.

Example 2 β€” fixing the rate that produced the inf

def safe_rate(zones, count_col, denom_col, *, min_denom=1, per=1000):
    """Compute a rate, making suppressed and zero denominators explicit."""
    counts = zones[count_col].to_numpy(dtype="float64")
    denom = zones[denom_col].to_numpy(dtype="float64")

    # official statistics often use negative sentinels for suppressed values
    suppressed = (counts < 0) | (denom < 0)
    too_small = denom < min_denom

    rate = np.full(len(zones), np.nan)
    usable = ~suppressed & ~too_small
    rate[usable] = counts[usable] / denom[usable] * per

    print(f"{usable.sum()} usable, {suppressed.sum()} suppressed, "
          f"{too_small.sum()} with denominator < {min_denom}")
    return rate, usable


zones["rate"], usable = safe_rate(zones, "count", "population", min_denom=50)
analysis = zones[usable].copy()
print(f"analysing {len(analysis)} of {len(zones)} zones")
63 usable, 2 suppressed, 1 with denominator < 50
analysing 63 of 66 zones

min_denom is doing real work beyond avoiding the division by zero. A zone with three residents and one event has a rate of 333 per 1,000, which will be the maximum on your map and is entirely an artefact of the denominator. Excluding tiny denominators is the crude fix; Empirical Bayes smoothing via esda.moran.Moran_Rate is the proper one.

Example 3 β€” comparing weights schemes when islands exist

import pandas as pd


def islands_sensitivity(zones, column):
    rows = []
    for label, scheme, kwargs in [
        ("queen, islands kept", "queen", {"islands": "report"}),
        ("queen, islands dropped", "queen", {"islands": "drop"}),
        ("knn k=4", "knn", {"k": 4}),
        ("knn k=8", "knn", {"k": 8}),
    ]:
        w, used = build_weights(zones, scheme, **kwargs)
        moran = safe_moran(used, column, w)
        rows.append({
            "weights": label,
            "n_used": w.n - len(w.islands),
            "I": round(moran.I, 4),
            "p": moran.p_sim,
        })
    frame = pd.DataFrame(rows)
    print(f"\nI spread: {frame['I'].max() - frame['I'].min():.4f}")
    return frame


print(islands_sensitivity(analysis, "rate").to_string(index=False))
                weights  n_used       I     p
    queen, islands kept      61  0.7877 0.001
 queen, islands dropped      61  0.7934 0.001
                knn k=4      63  0.7846 0.001
                knn k=8      63  0.7803 0.001

I spread: 0.0131

A spread of 0.013 across four treatments means the islands are not driving the result. If the spread were 0.3, the choice would be the finding, and it would need justifying in the text rather than buried in a default.

Explanation

Why nan comes with a plausible p-value

Moran computes the statistic and the permutation p-value through separate paths. The observed I is a ratio of sums over the values; one nan makes every sum nan.

The permutation p-value is the proportion of shuffled arrangements whose I was at least as extreme as the observed one. Comparisons against nan are always False in NumPy, so no permutation ever "beats" the observed value, and the count of exceedances stays at zero β€” which produces the smallest possible p-value, 1/(permutations+1) = 0.001.

So a broken statistic yields the most significant-looking p-value available. This is why np.isfinite(moran.I) belongs in the code and not in a comment.

Why islands are silently excluded rather than an error

Zero neighbours means zero weight in every term. Mathematically the unit simply drops out; there is no division by zero and nothing to raise about.

libpysal warns because the situation is usually unintended, but the warning is easy to lose in a notebook and impossible to see in a scheduled job that discards stderr. The number to report is w.n - len(w.islands) β€” the units that actually contributed.

The deeper issue is what an island means for your model. A genuine island has no land neighbours, so a contiguity model says nothing influences it. If your process is economic or social rather than physical, that is probably wrong, and KNN or a distance band is a better model of influence than shared boundaries.

Why disconnected components matter beyond islands

2 islands of 66 units
3 disconnected components

Three components, two of which are the single islands β€” so the mainland is one connected block. But a study area split into two large disconnected halves (an archipelago, a region separated by an international boundary you excluded) is a different situation: the statistic then measures autocorrelation within each half with no comparison between them.

Check w.n_components. More than one plus the island count means you have genuinely separate sub-regions, and a single global statistic across them is answering an odd question.

A single NaN value propagating through the sums to make Moran's I nan, while the permutation comparison yields the smallest possible p-value.
Comparisons against nan are always false, so no permutation beats the observed value and p comes out at its minimum.

Why "not fully connected" is sometimes a geometry bug

Not every island is a real island. Two polygons that should share a boundary may fail to touch because of:

  • precision β€” coordinates that differ in the eighth decimal place
  • slivers β€” a hairline gap between neighbours from different sources
  • invalid geometry β€” a self-intersection that breaks the touch test

If a unit you know is inland comes back as an island, check the geometry before changing the weights scheme:

suspect = zones.loc[w.islands]
print(suspect.geometry.is_valid)
neighbours = zones[zones.geometry.buffer(1).intersects(suspect.geometry.iloc[0])]
print(f"{len(neighbours) - 1} units within 1 m but not touching")

A unit with neighbours 1 m away but no shared boundary is a topology problem, not a spatial-weights problem, and the fix is to snap or clean the geometry.

Edge cases or notes

  • w.islands is a list of index labels, not positions. Use .drop(index=...), not .drop(...) with integers, unless your index is positional.
  • use_index=True is required by current libpysal when building from a GeoDataFrame, and keeps the weights aligned with your index.
  • KNN weights are asymmetric. A being B's nearest neighbour does not make B A's. That is fine for Moran's I, which symmetrises internally, but surprising if you inspect the matrix.
  • Row-standardising weights with islands leaves those rows summing to zero, not one β€” another reason they contribute nothing.
  • Moran_Rate applies an Empirical Bayes correction and is the right tool for rates with varying denominators.
  • Zero variance also breaks Geary's C and Getis-Ord. The check is worth applying once, before any spatial statistic.
  • A -99999 sentinel is finite, so np.isfinite will not catch it. Check the minimum of every numeric column from an official source.
  • w.n_components > 1 + len(w.islands) means genuinely separate sub-regions, which is worth handling explicitly.

FAQ

Does an island cause Moran's I to return nan?

No. Islands get zero weight and are silently excluded. nan comes from a non-finite value in the data or from zero variance.

Why is the p-value significant when I is nan?

Comparisons against nan are always false, so no permutation exceeds the observed value and the p-value comes out at its minimum, 1/(permutations+1). Always check np.isfinite(moran.I).

Should I drop islands or use KNN weights?

Either, deliberately. Drop them when a unit with no land neighbours genuinely has no spatial relationship; use KNN when influence travels by distance rather than shared boundaries. Report which you chose.

What causes a inf in a rate column?

A zero denominator. Guard with a minimum population, and remember that very small denominators produce extreme rates that are artefacts rather than findings.

A unit I know is inland is reported as an island. Why?

Its geometry does not quite touch its neighbours β€” precision, a sliver, or an invalidity. Fix the geometry rather than the weights.

What does "3 disconnected components" mean?

The adjacency graph is in three separate pieces. If that is more than the number of islands, you have genuinely separate sub-regions, and one global statistic across them may not be meaningful.

How do I report n honestly?

As w.n - len(w.islands) β€” the units that actually contributed to the statistic, not the number of rows in the frame.