Section 10.3 Local Getis-Ord
The maps above and the ones we introduced in Chapter 7 (produced through kernel density estimation) are helpful for visualising the varying intensity of crime across the study region of interest. However, as we noted above, it may be hard to discern random variation from clustering of crime incidents simply using these techniques. There are a number of statistical tools that have been developed to extract the signal from the noise, to determine what is "hot" from random variation. Many of these techniques were originally designed to detect clusters of disease.
In crime analysis it is common to use a technique developed by Arthur Getis and JK Ord in the early 1990s (see [Getis and Ord (1992)] and [Ord and Getis (1995)]), partly because it was implemented in the widely successful ArcGIS software tool as part of their hot spot analysis tool. This statistic provides a measure of spatial dependence at the local level. They are not testing homogeneity, they are testing dependence of a given attribute around neighbouring areas (which could be operationalised as micro-cells in a grid, as above, or could be some form of administrative areas, like wards, or census boundaries like LSOAs or census blocks).
Whereas the tests we covered in the previous chapter measure dependence at the global level (whereas areas are surrounded by alike areas), these measures of local clustering try to identify pockets of dependence. They are also called local indicators of spatial autocorrelation. It is indeed useful for us to be able to assess quantitatively whether crime events cluster in a non-random manner. In the words of Jerry [Ratcliffe (2010)], however, "while a global Moran’s I test can show that crime events cluster in a non-random manner, this simply explains what most criminal justice students learn in their earliest classes." (p. 13). For a crime analyst and practitioner, while clustering is of interest, learning about the existence and location of local clusters is of paramount importance.
There are two variants of the local Getis-Ord statistic. The local \(G\) is computed as a ratio of the weighted average of the values of the attribute of interest in the neighbouring locations, not including the value at the location. It is generally used for spread and diffusion studies. The local \(G^*\) on the other hand, includes the value at the location in both numerator and denominator. It is more generally used for clustering studies (for formula details, see [Getis and Ord, 1992] and [Ord and Getis, 1995]). High positive values indicate the possibility of a local cluster of high values of the variable being analysed ("hot spot"); very low relative values indicate a similar cluster of low values ("cold spot").
A value larger than the mean (or, a positive value for a standardized z-value) suggests a High-High cluster or hot spot, a value smaller than the mean (or, negative for a z-value) indicates a Low-Low cluster or cold spot ([Anselin, 2020]).
For statistical inference, a Bonferroni-type test is suggested in the papers by Getis and Ord, where tables of critical values are included. The critical values for the 95th percentile under the assumptions discussed by [Ord and Getis (1995)] are for \(n\)=30: 2.93, \(n\)=50: 3.08, \(n\)=100: 3.29, \(n\)=500: 3.71, and for \(n\)=1000: 3.89. [Anselin (2020)] suggests that this analytical approximation may not be reliable in practice and that conditional permutation is advisable. Inference with these local measures are also complicated by problems of multiple comparisons, so it is typically advisable to use some form of adjustment to this problem.
The first step in computing the Getis and Ord statistics involves defining the neighbours of each area. And to do that, first we have to turn each cell in our grid into a polygon, so we basically move from the raster to the vector model. This way we can use the methods discussed in Chapter 9 to operationalise neighbouringness in our data. To transform our raster grid cells into polygons, we can use the
rasterToPolygons() function from the raster package.
# Move to vector model
grid_cc <- rasterToPolygons(burglary_raster_cc)
What is a neighbour and the code related to defining neighbours is a topic we covered in Chapter 9. For brevity, we will only use the queen criteria for this illustration. Only to indicate that
include.self() ensures we get the \(G^*\) variation of the local Getis-Ord test, the second variation we mentioned above.
nb_queen_cc <- poly2nb(grid_cc)
lw_queen_cc <- nb2listw(include.self(nb_queen_cc), zero.policy=T)
Now that we have the data ready we can, analytically, compute the statistic using the
localG() function from the spdep package. We can then add the computed local \(Gi^*\) to each cell. And we can then see the produced z scores:
# Perform the local G analysis (Getis-Ord GI*)
grid_cc$G_hotspot_z <- as.vector(localG(grid_cc$layer, lw_queen_cc))
# Summarise results
summary(grid_cc$G_hotspot_z)
## Min. 1st Qu. Median Mean 3rd Qu. Max. ## -1.411 -1.011 -0.398 0.033 0.835 5.231
How do you interpret these? They are z scores. So you expect large absolute values to index the existence of cold spots (if they are negative) or hot spots (if they are positive). Our sample size is 475 cells, and there are only a very small handful of cells with values that reach the critical value in our dataset.
A known problem with this statistic is that of multiple comparisons. As noted by [Pebesma and Bivand (2021)]: "although the apparent detection of hot spots from values of local indicators has been quite widely adopted, it remains fraught with difficulty because adjustment of the inferential basis to accommodate multiple comparisons is not often chosen." So we need to ensure the critical values we use to adjust for this. As recommended by [Pebesma and Bivand (2021)], we can use the
p.adjust() function to count the number of observations that this statistic identifies as hot spots using different solutions for the multiple comparison problem, such as Bonferroni, False Discovery Rate (fdr), or the method developed by Benjamini and Yekutielli (BY).
p_values_lg <- 2 * pnorm(abs(c(grid_cc$G_hotspot_z)), lower.tail = FALSE)
p_values_lgmc <- cbind(p_values_lg,
p.adjust(p_values_lg, "bonferroni"),
p.adjust(p_values_lg, "fdr"),
p.adjust(p_values_lg, "BY"))
colnames(p_values_lgmc) <- c("None", "Bonferroni",
"False Discovery Rate", "BY (2001)")
apply(p_values_lgmc, 2, function(x) sum(x < 0.05))
## None Bonferroni ## 53 10 ## False Discovery Rate BY (2001) ## 15 12
We can see that we get fewer significant hot (or cold) spots once we adjust for the problem of multiple comparisons. In general, Bonferroni correction is the more conservative option. Here, the p-value is multiplied by the number of comparisons. This is useful, when let’s say we’re trying to avoid false positives. We will progress with those areas identified with this method.
To map this, we can create a new variable, let’s call it
hotspot, which can take the value of "hot spot" if the Bonferroni adjusted p-value is smaller than 0.05 and the z score is positive, and "cold spot" if the Bonferroni adjusted p-value is smaller than 0.05 and the z score is negative. Otherwise, it will be labelled "Not significant".
# convert to simple features object
grid_cc_sf <- st_as_sf(grid_cc)
# Creates p values for Bonferroni adjustment
grid_cc_sf$G_hotspot_bonf <- p.adjust(p_values_lg, "bonferroni")
# Create character vector based on Z score and p values
grid_cc_sf <- grid_cc_sf %>%
mutate(hotspot = case_when(grid_cc_sf$G_hotspot_bonf < 0.05 &
grid_cc_sf$G_hotspot_z > 0 ~ "Hot spot",
grid_cc_sf$G_hotspot_bonf < 0.05 &
grid_cc_sf$G_hotspot_z < 0 ~ "Cold spot",
TRUE ~ "Not significant"))
Now let’s see the result:
table(grid_cc_sf$hotspot)
## ## Hot spot Not significant ## 10 512
Great, we have our 10 significant hot spots, and no significant cold spots apparently. We can now map these values.
#set to view mode and make sure colours are colourblind safe
tmap_mode("view")
tmap_style("col_blind")
# make map of the z-scores
map1 <- tm_shape(grid_cc_sf) +
tm_fill(c("G_hotspot_z"), title="Z scores", alpha = 0.5)
# make map of significant hot spots
map2 <- tm_shape(grid_cc_sf) +
tm_fill(c("hotspot"), title="Hot spots", alpha = 0.5)
# print both maps side by side
tmap_arrange(map1, map2, nrow=1)
We can see here that not all the micro-grid cells with an elevated count of burglaries exhibit significant local dependence. There are only two clusters of burglaries where we see positive spatial dependence: one is located north of Picadilly Gardens in the Northern Quarter of Manchester (an area of trendy bars and shops, and significant gentrification, in proximity to a methadone delivery clinic) and the Deansgate area of the city centre (one of the main commercial hubs of the city not far from the criminal courts).
