Scribd
Upload
29 views10 pages

Spatial Autocorrelation Analysis in R

This document provides instructions for running measures of spatial autocorrelation in R. It will cover calculating global and local spatial autocorrelation statistics. First, necessary libraries are loaded and shapefiles and data are imported. Then instructions are given to calculate neighbors and create a weights list. A global Moran's I test is run to measure overall spatial autocorrelation. Local indicators of spatial autocorrelation are also calculated and mapped to identify clustering in different areas.

Uploaded by

Med Abdel
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd
29 views10 pages

Spatial Autocorrelation Analysis in R

This document provides instructions for running measures of spatial autocorrelation in R. It will cover calculating global and local spatial autocorrelation statistics. First, necessary libraries are loaded and shapefiles and data are imported. Then instructions are given to calculate neighbors and create a weights list. A global Moran's I test is run to measure overall spatial autocorrelation. Local indicators of spatial autocorrelation are also calculated and mapped to identify clustering in different areas.

Uploaded by

Med Abdel
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd

Practical 9: Measuring Spatial

Autocorrelation in R
An Introduction to Spatial Data Analysis and Visualisation in R - Guy Lansley & James Cheshire (2016)

This practical will cover how to run various measures of spatial autocorrelation in R. We will consider both
statistics of global spatial autocorrelation and how to identify spatial clustering across our study area. Data for
the practical can be downloaded from the Introduction to Spatial Data Analysis and Visualisation in R
([Link] homepage.

In this practical we will:

Run a global Spatial autocorrelation for a shapefile


Identify local indicators of spatial autocorrelation
Run a Getis-Ord

First, we must set the working directory and load the practical data.

# Set the working directory


setwd("C:/Users/Guy/Documents/Teaching/CDRC/Practicals")

# Load the data. You may need to alter the file directory
[Link] <-[Link]("practical_data.csv")

We will also need to load the spatial data files from the previous practicals.

# load the spatial libraries


library("sp")
library("rgdal")
library("rgeos")

# Load the output area shapefiles


[Link] <- readOGR(".", "Camden_oa11")

## OGR data source with driver: ESRI Shapefile


## Source: ".", layer: "Camden_oa11"
## with 749 features
## It has 1 fields

# join our census data to the shapefile


[Link] <- merge([Link], [Link], by.x="OA11CD", by.y="OA")

# load the houses point files


[Link] <- readOGR(".", "Camden_house_sales")

## OGR data source with driver: ESRI Shapefile


## Source: ".", layer: "Camden_house_sales"
## with 2547 features
## It has 4 fields

Remember the distribution of our qualification variable? We will be working on that today. We have first
mapped it to remind us of its spatial distribution across our study area.
library("tmap")

tm_shape([Link]) + tm_fill("Qualification", palette = "Reds", style = "quantile", title =


"% with a Qualification") + tm_borders(alpha=.4)

Running a spatial autocorrelation


A spatial autocorrelation ([Link] measures
how distance influences a particular variable. In other words, it quantifies the degree of which objects are
similar to nearby objects. Variables are said to have a positive spatial autocorrelation when similar values tend
to be nearer together than dissimilar values.

Waldo Tober’s first law of geography is that “Everything is related to everything else, but near things are more
related than distant things.” so we would expect most geographic phenomena to exert a spatial autocorrelation
of some kind. In population data this is often the case as persons with similar characteristics tend to reside in
similar neighbourhoods due to a range of reasons including house prices, proximity to workplaces and cultural
factors.

We will be using the spatial autocorrelation functions available from the spdep package.

library(spdep)

Finding neighbours
In order for the subsequent model to work, we need to work out what polygons neighbour each other. The
following code will calculate neighbours for our [Link] polygon and print out the results below.
# Calculate neighbours
neighbours <- poly2nb([Link])
neighbours

## Neighbour list object:


## Number of regions: 749
## Number of nonzero links: 4342
## Percentage nonzero weights: 0.7739737
## Average number of links: 5.797063

We can plot the links between neighbours to visualise their distribution across space.

plot([Link], border = 'lightgrey')


plot(neighbours, coordinates([Link]), add=TRUE, col='red')

# Calculate the Rook's case neighbours


neighbours2 <- poly2nb([Link], queen = FALSE)
neighbours2

## Neighbour list object:


## Number of regions: 749
## Number of nonzero links: 4176
## Percentage nonzero weights: 0.7443837
## Average number of links: 5.575434

We can already see that this approach has identified fewer links between neighbours. By plotting both
neighbour outputs we can interpret their differences.
# compares different types of neighbours
plot([Link], border = 'lightgrey')
plot(neighbours, coordinates([Link]), add=TRUE, col='blue')
plot(neighbours2, coordinates([Link]), add=TRUE, col='red')

We can represent spatial autocorrelation in two ways; globally or locally. Global models
([Link] will create a single measure
which represents the entire data whilst local models ([Link]
local_indicators_of_spatial_as.htm) let us explore spatial clustering across space.

Running a global spatial autocorrelation


With the neighbours defined. We can now run a model. First, we need to convert the data types of the
neighbours object. This file will be used to determine how the neighbours are weighted

# Convert the neighbour data to a listw object


listw <- nb2listw(neighbours2)
listw
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 749
## Number of nonzero links: 4176
## Percentage nonzero weights: 0.7443837
## Average number of links: 5.575434
##
## Weights style: W
## Weights constants summary:
## n nn S0 S1 S2
## W 749 561001 749 285.3793 3113.982

We can now run the model. This type of model is known as a Moran’s test. This will create a correlation score
between -1 and 1. Much like a correlation coefficient, 1 determines perfect positive spatial autocorrelation (so
our data is clustered), 0 identifies the data is randomly distributed and -1 represents negative spatial
autocorrelation (so dissimilar values are next to each other).

# global spatial autocorrelation


[Link]([Link]$Qualification, listw)

##
## Moran I test under randomisation
##
## data: [Link]$Qualification
## weights: listw
##
## Moran I statistic standard deviate = 24.292, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.5448699398 -0.0013368984 0.0005055733

The Moran I statistic is 0.54, we can, therefore, determine that there our qualification variable is positively
autocorrelated in Camden. In other words, the data does spatially cluster. We can also consider the p-value as
a measure of the statistical significance of the model.

Running a local spatial autocorrelation


We will first create a moran plot which looks at each of the values plotted against their spatially lagged values.
It basically explores the relationship between the data and their neighbours as a scatter plot. The style refers to
how the weights are coded. “W” weights are row standardised (sums over all links to n).

# creates a moran plot


moran <- [Link]([Link]$Qualification, listw = nb2listw(neighbours2, style = "W"))
Is it possible to determine a positive relationship from observing the scatter plot?

# creates a local moran output


local <- localmoran(x = [Link]$Qualification, listw = nb2listw(neighbours2, style = "W"))

By considering the help page for the localmoran function (run ?localmoran in R) we can observe the
arguments and outputs. We get a number of useful statistics from the model which are as defined:

Name Description

Ii local moran statistic

[Link] expectation of local moran statistic

[Link] variance of local moran statistic

[Link] standard deviate of local moran statistic

Pr() p-value of local moran statistic

First, we will map the local moran statistic (Ii). A positive value for Ii indicates that the unit is surrounded by
units with similar values.

# binds results to our polygon shapefile


[Link] <- cbind([Link], local)

# maps the results


tm_shape([Link]) + tm_fill(col = "Ii", style = "quantile", title = "local moran statisti
c")
From the map, it is possible to observe the variations in autocorrelation across space. We can interpret that
there seems to be a geographic pattern to the autocorrelation. However, it is not possible to understand if
these are clusters of high or low values.

Why not try to make a map of the P-value to observe variances in significance across Camden? Use
names([Link]@data) to find the column headers.

One thing we could try to do is to create a map which labels the features based on the types of relationships
they share with their neighbours (i.e. high and high, low and low, insignificant, etc…). The following code will
run this for you. Source: Brunsdon and Comber (2015) ([Link]
for-spatial-analysis-and-mapping/book241031)
### to create LISA cluster map ###
quadrant <- vector(mode="numeric",length=nrow(local))

# centers the variable of interest around its mean


[Link] <- [Link]$Qualification - mean([Link]$Qualification)

# centers the local Moran's around the mean


[Link] <- local[,1] - mean(local[,1])

# significance threshold
signif <- 0.1

# builds a data quadrant


quadrant[[Link] >0 & [Link]>0] <- 4
quadrant[[Link] <0 & [Link]<0] <- 1
quadrant[[Link] <0 & [Link]>0] <- 2
quadrant[[Link] >0 & [Link]<0] <- 3
quadrant[local[,5]>signif] <- 0

# plot in r
brks <- c(0,1,2,3,4)
colors <- c("white","blue",rgb(0,0,1,alpha=0.4),rgb(1,0,0,alpha=0.4),"red")
plot([Link],border="lightgray",col=colors[findInterval(quadrant,brks,[Link]=FALSE)])
box()
legend("bottomleft",legend=c("insignificant","low-low","low-high","high-low","high-high"),
fill=colors,bty="n")

It is apparent that there is a statistically significant geographic pattern to the clustering of our qualification
variable in Camden.
Getis-Ord
Another approach we can take is hot-spot analysis. The Getis-Ord Gi Statistic looks at neighbours within a
defined proximity to identify where either high or low values cluster spatially. Here statistically significant hot-
spots are recognised as areas of high values where other areas within a neighbourhood range also share high
values too.

First, we need to define a new set of neighbours. Whilst the spatial autocorrection considered units which
shared borders, for Getis-Ord we are defining neighbours based on proximity. The example below shows the
results where we have a search radius of 250 metres.

However, here a search radius of just 250 metres fails to define nearest neighbours for some areas so we will
need to set the radius as 800 metres or more for our model in Camden.

# creates centroid and joins neighbours within 0 and 800 units


nb <- dnearneigh(coordinates([Link]),0,800)
# creates listw
nb_lw <- nb2listw(nb, style = 'B')

# plot the data and neighbours


plot([Link], border = 'lightgrey')
plot(nb, coordinates([Link]), add=TRUE, col = 'red')

With a set of neighbourhoods established we can now run the test and bind the results to our polygon file.

On some machines the cbind may not work with a spatial data file, in this case, you will need to change
[Link] to [Link]@data ([Link] so that R knows which part of the spatial data file
to join. If you take this approach the subsequent column ordering may be different to what is shown in the
example below.

# compute Getis-Ord Gi statistic


local_g <- localG([Link]$Qualification, nb_lw)
local_g <- cbind([Link], [Link](local_g))
names(local_g)[6] <- "gstat"

# map the results


tm_shape(local_g) + tm_fill("gstat", palette = "RdBu", style = "pretty") + tm_borders(alpha=.
4)

The Gi Statistic is represented as a Z-score. Greater values represent a greater intensity of clustering and the
direction (positive or negative) indicates high or low clusters. The final map should indicate the location of hot-
spots across Camden. Repeat this for another variable.

The rest of the online tutorials in this series can be found at: [Link]
spatial-data-analysis-and-visualisation-r ([Link]
visualisation-r)

Common questions

Powered by AI

Spatial autocorrelation is significant in geographical studies because it quantifies how much a variable is similar between nearby spatial locations. This helps in understanding patterns and clusters within geographical data, such as why similar demographic characteristics might be located near each other due to factors like house prices, proximity to workplaces, and cultural elements . It also aids in identifying hotspots or areas with significant clustering, which is essential for urban planning and resource allocation .

Using both global and local models is essential when analyzing spatial autocorrelation to gain a comprehensive understanding of spatial patterns. Global models like Moran's I provide an overall measure indicating whether clustering exists across the entire study region, while local models like LISA provide detailed insight into where and how these patterns occur at a local level. This dual approach ensures that subtle and localized spatial patterns are not overlooked .

Handling spatial data in R, specifically polygon files, can present challenges such as ensuring accurate projection and alignment of spatial layers, dealing with large datasets that strain system resources, and requiring correct manipulation of data types for various functions. Additionally, operations like merging data and calculating spatial relationships require precise coding, and even minor errors can lead to incorrect interpretations or software errors, such as failing to properly cbind spatial data to non-spatial data .

The Getis-Ord statistic distinguishes itself from other methods of spatial analysis by focusing on identifying clusters of high or low values within a defined proximity, known as hotspots or cold spots. This is different from methods like Moran's I which measure overall autocorrelation. The Getis-Ord statistic uses a Z-score to represent intensity and direction of clustering, enabling the identification of specific areas that significantly deviate from the norm .

A lower percentage of nonzero weights in the neighbor list implies fewer spatial connections between geographic units, which can reduce the sensitivity of detecting spatial patterns and clustering. It could lead to higher rates of false negatives in identifying significant relationships. In spatial studies, maintaining an optimal percentage of nonzero weights is critical for ensuring robustness and reliability in resulting analyses, such as when fewer links are observed in the Rook’s case compared to the Queen’s neighbor definition .

It is necessary to plot both Pearson and Rook’s case neighbors to compare their definitions of proximity and assess their impact on spatial analysis. Pearson adjacency considers sharing a common side (edge) as defining neighbors, whereas Rook’s case considers sharing a vertex or border. This comparison helps in understanding how different neighborhood conceptualizations influence results, such as in spatial clustering and pattern detection, as they offer varying sensitivity to spatial configurations .

The Moran's I statistic measures spatial autocorrelation by assessing the correlation of a variable with itself across a spatial area. A Moran's I value of 1 indicates a perfect positive spatial autocorrelation (clusters of similar values), 0 denotes random distribution, and -1 indicates perfect negative autocorrelation (dissimilar values adjacent to each other). In the context of Camden's qualification data, a Moran's I of 0.54 suggests a positive spatial autocorrelation, meaning the data is spatially clustered .

LISA, or Local Indicators of Spatial Association, plays a crucial role in spatial data analysis by allowing the exploration of spatial clusters at a local level, rather than providing a single summary measure as global models do. It enables the identification of spatial patterns, showing where clusters of high or low values occur, and their statistical significance, thus providing insights into local spatial variations and dependencies in the data .

Choosing an appropriate distance is crucial in the Getis-Ord Gi statistic to accurately identify spatial clusters. If the distance is too small, it may fail to capture relevant neighbouring areas, leading to an incomplete analysis. Conversely, a distance that is too large may include non-relevant areas, diluting significant clusters. For example, in analyzing Camden's data, a radius of 250 meters was initially too small, necessitating an increase to 800 meters for effective identification of spatial hot and cold spots .

The choice of weighting style, such as row-standardized weights ("W" weights), impacts spatial autocorrelation analysis by determining how relationships between spatial units are modeled in the analysis. "W" weights standardize the influence of each neighbor, which affects the calculation of statistics like Moran's I and can influence the detection of spatial patterns. A different weighting style could potentially yield different insights, affecting interpretations of cluster significance in spatial datasets .

You might also like