Analyze Spatial Clusters and Hot Spots Using PySAL and Python
Spatial data analysis is a crucial part of Geographic Information Systems (GIS), remote sensing, urban planning, environmental observation, and location intelligence. Instead of being scattered randomly, geographical phenomena often have particular spatial patterns. The identification of such patterns enables scholars as well as professionals in GIS to learn about the relationships between geographical objects and areas with high concentrations of population.
PySAL (Python Spatial Analysis Library) is an open-source Python library specifically developed for spatial data science, spatial statistics, and geospatial modeling. It contains tools to analyze spatial autocorrelation, recognize spatial clusters, identify significant hot and cold spots, and analyze geographical correlation among different spatial datasets.
By means of PySAL and with the help of Python libraries such as GeoPandas, NumPy, Matplotlib, and SciPy, GIS analysts can perform advanced statistical analysis of spatial data and find geographical patterns.
What Is Spatial Cluster Analysis in GIS?
Spatial cluster analysis is a type of statistical analysis that is utilized for the identification of locations where some traits or values are clustered. Spatial clustering is different from the standard clustering procedures in that it also takes into account the geographical relations between the studied objects.
For example, spatial cluster analysis is necessary to identify the following:
Areas with high population density;
Places with disease outbreaks;
The spatial distribution of pollution;
Urban areas with high crime rates;
Agricultural areas with similar types of vegetation;
Areas where land usage changes.
Spatial clustering is based on the principle of spatial autocorrelation that assists in identifying whether the given observations are more or less similar in relation to the distance between them.
Positive spatial autocorrelation indicates that two neighboring areas are likely to have the same trait. Negative spatial autocorrelation states that two neighboring areas show different traits.
What Is Hot Spot Analysis?
Hot spot analysis is a way to locate geographic areas where there are noteworthy amounts of certain qualities and attributes higher than expected based on the analysis done with the help of statistical methods for spatial distribution.
In order to identify hot or cold spots, the hot spot analysis defines the areas where certain qualities have statistically significant spatial concentration.
For example, in a study of urban heat islands, hot spots will be where neighborhoods have significantly higher surface temperatures than those elsewhere, while cold spots will be areas with colder temperatures than average.
The hot spot analysis is commonly utilized in different areas of GIS applications, such as environmental monitoring, urban planning, public health, and so on.
What is PySAL in Python?
PySAL is a set of Python packages for spatial analysis, geospatial statistics, and spatial modeling.
PySAL provides a modular framework that enables GIS specialists to carry out geographic data analysis using statistical techniques.
Some of its major components include:
libpysal: Provides spatial weights, computational geometry, and supporting data structures.
esda: Offers exploratory spatial data analysis statistics, including Moran's I, local Moran's I, and Getis-Ord statistics.
spreg: Offers spatial regression models for geographic relationship analysis.
giddy: Offers support for spatial dynamics and changes in geographic distribution over time.
mapclassify: Provides thematic mapping classification.
How to Install PySAL and Required Python Libraries
Before conducting spatial statistical analysis, install the required Python packages.
pip install pysal geopandas matplotlib numpy scipyAlternatively, install only the primary packages needed for this workflow:
pip install libpysal esda geopandas matplotlibImport the necessary modules:
import geopandas as gpd
import matplotlib.pyplot as plt
import numpy as np
from libpysal.weights import Queen
from esda.moran import Moran, Moran_Local
from esda.getisord import G_Local
These libraries support geospatial data manipulation, spatial weights construction, spatial autocorrelation calculations, and geographic visualization.
Step 1: Load Geospatial Data Using GeoPandas
GeoPandas provides a convenient interface for reading and manipulating geographic datasets in Python.
It supports multiple geospatial formats, including Shapefile, GeoJSON, and GeoPackage.
import geopandas as gpd
gdf = gpd.read_file("study_area.geojson")
print(gdf.head())
print(gdf.crs)
print(gdf.columns)
Assume the dataset contains administrative boundaries and an attribute called population_density.
gdf = gdf[
gdf.geometry.notna()
& gdf["population_density"].notna()
].copy()
gdf = gdf.reset_index(drop=True)
This removes records with missing geometries or attribute values.
For reliable spatial analysis, geometries should also be checked for validity, and geographic units should be appropriate for the analysis.
Step 2: Create Spatial Weights Using PySAL
Spatial weights describe how geographic features are related to their neighbors.
Spatial weights matrices are fundamental to spatial statistical analysis because they define the geographic relationships used to calculate spatial autocorrelation.
Common spatial weights methods include:
Queen contiguity: Polygons sharing a boundary or vertex are considered neighbors.
Rook contiguity: Polygons sharing a boundary segment are considered neighbors.
Distance-based weights: Features within a specified distance are considered neighbors.
K-nearest neighbors: Each feature is connected to a specified number of nearby features.
Create Queen Contiguity Spatial Weights
from libpysal.weights import Queen
w = Queen.from_dataframe(
gdf,
use_index=False
)
w.transform = "r"
print("Number of observations:", w.n)
print("Number of neighbor links:", w.s0)
The r transformation row-standardizes the spatial weights matrix so that weights for observations with neighbors sum to one.
Queen contiguity is commonly used for polygon datasets such as census tracts, administrative districts, and land parcels.
Important: Check for islands, or polygons without neighbors. These may require a different neighborhood definition, such as distance-based or nearest-neighbor weights.
Step 3: Calculate Global Moran's I
Global Moran's I measures the overall degree of spatial autocorrelation in a geographic dataset.
It determines whether high or low values tend to occur near similar values across the entire study area.
The mathematical formula is:
\[ I=\frac{n}{S_0}\frac{\sum_i\sum_j w_{ij}(x_i-\bar{x})(x_j-\bar{x})}{\sum_i(x_i-\bar{x})^2} \]
Where:
\(n\) = number of spatial observations.
\(w_{ij}\) = spatial weight between locations \(i\) and \(j\).
\(x_i\) = attribute value at location \(i\).
\(\bar{x}\) = mean attribute value.
\(S_0\) = sum of all spatial weights.
Python Implementation
from esda.moran import Moran
y = gdf["population_density"].to_numpy()
moran = Moran(
y,
w,
permutations=999
)
print("Moran's I:", moran.I)
print("Permutation p-value:", moran.p_sim)
print("Z-score:", moran.z_sim)
How to Interpret Moran's I
A positive Moran's I indicates positive spatial autocorrelation, meaning neighboring geographic features tend to have similar values.
A negative Moran's I indicates spatial dispersion or neighboring dissimilarity.
Values near the statistic's expected value under the null hypothesis suggest little global spatial autocorrelation.
The permutation-based p-value helps determine whether the observed spatial pattern is statistically significant. Moran's I should not be interpreted solely by comparing its value with zero because its expected value depends on the number of observations and the spatial weights structure.
Step 4: Identify Local Spatial Clusters Using Local Moran's I
Global Moran's I measures overall spatial autocorrelation but does not identify the specific geographic locations contributing to spatial clustering.
Local Moran's I, also known as a Local Indicator of Spatial Association (LISA), helps detect local clusters and spatial outliers.
The method commonly identifies four spatial association categories:
Category | Description |
High-High (HH) | High-value feature surrounded by high-value neighbors |
Low-Low (LL) | Low-value feature surrounded by low-value neighbors |
High-Low (HL) | High-value feature surrounded by low-value neighbors |
Low-High (LH) | Low-value feature surrounded by high-value neighbors |
Calculate Local Moran's I
from esda.moran import Moran_Local
local_moran = Moran_Local(
y,
w,
permutations=999,
seed=42
)
gdf["local_I"] = local_moran.Is
gdf["p_value"] = local_moran.p_sim
gdf["quadrant"] = local_moran.q
The q attribute classifies observations into spatial association quadrants.
Using PySAL's default quadrant coding:
1 = High-High
2 = Low-High
3 = Low-Low
4 = High-Low
Classify Significant Spatial Clusters
labels = {
1: "High-High",
2: "Low-High",
3: "Low-Low",
4: "High-Low"
}
gdf["cluster"] = "Not Significant"
significant = gdf["p_value"] < 0.05
gdf.loc[significant, "cluster"] = (
gdf.loc[significant, "quadrant"]
.map(labels)
This creates a categorical attribute identifying statistically significant local spatial relationships.
A 0.05 significance threshold is commonly used, but results from many local tests should be interpreted carefully because multiple comparisons can increase false-positive findings.
Step 5: Visualize Spatial Clusters in Python
Visualizing spatial clusters helps analysts identify geographic areas exhibiting similar spatial patterns.
GeoPandas and Matplotlib can generate thematic maps from the Local Moran's I results.
import matplotlib.pyplot as plt
fig, ax = plt.subplots(
figsize=(12, 8)
)
gdf.plot(
column="cluster",
categorical=True,
legend=True,
ax=ax,
edgecolor="black",
linewidth=0.2
)
ax.set_title(
"Local Moran's I Spatial Clusters"
)
ax.set_axis_off()
plt.tight_layout()
plt.show()
The resulting map highlights statistically significant high-high and low-low clusters, along with high-low and low-high spatial outliers.
For example, high-high clusters of population density may indicate urban centers, whereas low-low clusters may indicate less densely populated geographic regions.
Step 6: Perform Hot Spot Analysis Using Getis-Ord Gi*
Getis-Ord Gi* is a local spatial statistic used to identify statistically significant concentrations of high and low attribute values.
Unlike Local Moran's I, which evaluates local similarity and spatial outliers, Getis-Ord Gi* focuses on concentrations of high or low values within defined neighborhoods.
The Gi* statistic is based on a weighted sum of attribute values around each location, including the location itself.
Python Implementation
For Gi* analysis, construct binary spatial weights and explicitly include each feature in its own neighborhood.
from libpysal.weights import Queen
from esda.getisord import G_Local
w_gi = Queen.from_dataframe(
gdf,
use_index=False
)
w_gi.transform = "b"
hotspots = G_Local(
y,
w_gi,
transform="B",
star=True,
permutations=999,
seed=42
)
gdf["gi_zscore"] = hotspots.Zs
gdf["gi_pvalue"] = hotspots.p_sim
Here, star=True requests Gi* rather than the ordinary Gi statistic. The binary transformation and explicit self-inclusion ensure that the neighborhood definition matches the intended Gi* analysis.
Identify Hot Spots and Cold Spots
gdf["hotspot_type"] = "Not Significant"
significant = gdf["gi_pvalue"] < 0.05
gdf.loc[
significant &
(gdf["gi_zscore"] > 0),
"hotspot_type"
] = "Hot Spot"
gdf.loc[
significant &
(gdf["gi_zscore"] < 0),
"hotspot_type"
] = "Cold Spot"
This code classifies significant local concentrations based on permutation p-values and the direction of the standardized Gi* statistic.
The z-score indicates the direction of local concentration, while the permutation-based p-value provides the significance assessment.
Step 7: Create a Hot Spot Map
fig, ax = plt.subplots(
figsize=(12, 8)
)
gdf.plot(
column="hotspot_type",
categorical=True,
legend=True,
ax=ax,
edgecolor="gray",
linewidth=0.2
)
ax.set_title(
"Getis-Ord Gi* Hot Spot Analysis"
)
ax.set_axis_off()
plt.tight_layout()
plt.show()
The output map identifies geographic areas containing statistically significant concentrations of high and low values.
Hot spot maps can be used to support decisions about resource allocation, environmental monitoring, infrastructure planning, and geographic risk assessment.
The use of spatial cluster detection and hot spot analysis is of huge importance in the analysis of geographical patterns in different data sets of GIS. PySAL serves as a very convenient tool for applying spatial statistics to geographical data since the PySAL library is based on the Python programming language.
GIS specialists can conduct measurements of spatial autocorrelation with the help of Global Moran's I, Local Moran's I, and Getis-Ord Gi* methods that enable them to learn about the existence of geographical clusters with high or low statistical significance and spatial outliers on the map.
Moreover, when PySAL is used together with GeoPandas and remote sensing datasets, as well as some tools for visualization, it becomes easy for developers and scientists to conduct their research on different spheres of scientific life, such as urban partitioning, environmental research, public health, agriculture, and location intelligence.
To learn more about PySAL and its geospatial capabilities, click here.
For more information or any questions regarding the LizardTech suite of products, please don't hesitate to contact us at:
Email: info@geowgs84.com
USA (HQ): (720) 702–4849
(A GeoWGS84 Corp Company)
Comments