1 Overview

1.1 Context

Since late December 2019, an outbreak of a novel coronavirus disease (COVID-19; previously known as 2019-nCoV) was reported in Wuhan, China, which has subsequently affected 210 countries worldwide. In general, COVID-19 is an acute resolved disease but it can also be deadly, with a 2% case fatality rate.

The COVID-19 pandemic in Mexico is part of the ongoing worldwide pandemic of coronavirus disease 2019 (COVID-19) caused by severe acute respiratory syndrome coronavirus 2 (SARSCoV-2). The virus was confirmed to have reached Mexico in February 2020. However, the National Council of Science and Technology (CONACYT) reported two cases of COVID-19 in midJanuary 2020 in the states of Nayarit and Tabasco, one case per state. As of September 26, there had been 726,431 confirmed cases of COVID-19 in Mexico and 76,243 reported deaths, although the Secretariat of Health, through the “Programa Centinela” (Spanish for “Sentinel Program”) estimated in mid July 2020 that there were more than 2,875,734 cases in Mexico, because they were considering the total number of cases confirmed as a statistical sample. (Note: This paragraph is extracted from wiki: https://en.wikipedia.org/wiki/COVID19_pandemic_in_Mexico).

1.2 Task and Objective

In this take-home exercise, the objective is to analyse the spatio-temporal patterns of COVID-19 case at Central Mexico (i.e. Mexico City (9), Mexico State (15) and Morelos State (17) by using localized spatial statistics methods.

1.3 Guide and Resources

The analysis techniques and methods will be follow along using the Prof’s Kam Hands On Exercise Materials as follow: - IS415 In-class Exercise 4 - https://rpubs.com/tskam/IS415-Hands-on_Ex06 - https://rpubs.com/tskam/IS415-Hands-on_Ex07

1.4 The Data Set

For the purpose of this study, a GIS data called municipalities_COVID is given. It is in ESRI shapefile format. The shapefile consists of reported COVID-19 cases and other related data at the municipality level. The data was extracted from Secretary of Health (Spanish: Secretaría de Salud), Government of Mexico homepage’s (http://datosabiertos.salud.gob.mx/gobmx/salud/datos_abiertos/).

The original data are stored in csv file format. The cases are recorded daily from 1st January 2020 until 6th August 2020. They had been aggregated up to epidemiological week (e-week) and municipality levels. There are a total of 32 weeks of records in the municipalities_COVID GIS data.

In this study, only required to focus on Mexico City, Mexico State and Morelos State. Their CVE_ENT codes are 9, 15 and 17 respectively. In terms of epidemiological period, and only focus on e-week 13 to e-week 32.

2 Setting Up Working Environment

This code chunk will check if the R packages in the packaging list have been installed. if no, go ahead to install the missing one. After the installation, it will also load the R packages in R.

#c is creating a list to a packages
#"<-" && "=" serves as the same function

packages <- c('sp', 'spdep', 'rgdal', 'rgeos', 'raster', 'sf', 'tmap', 'tidyverse', 'readr', 'ggplot2', 'spatstat', 'dplyr')
for (p in packages){
  if(!require(p, character.only = T)){
    install.packages(p)
  }
  library(p, character.only = T)
}

3 Geospatial Data Wrangling

3.1 Loading up the Geospatial Data

The code chunk below uses st_read() function of sf package to import the datas into R as simple feature dataframes. And also import data by using read_csv() of readr package. Read_excel is being used for the excel format of “xls” & “xlsx”. While “csv” uses read_csv.

municipalities_COVID = readOGR(dsn = "data/geospatial", layer="municipalities_COVID")
## OGR data source with driver: ESRI Shapefile 
## Source: "D:\Google Drive\SMU\Year 4\Term 1\IS415-Geospatial Analytics and Applications\Take Home Exercises\Take-home_Ex02\data\geospatial", layer: "municipalities_COVID"
## with 2465 features
## It has 198 fields
#municipalities_COVID_sf <- st_as_sf(municipalities_COVID)

# #Alternatively
#municipalities_COVID_sf = st_read(dsn = "data/geospatial", "municipalities_COVID")

#print(municipalities_COVID_sf)

3.2 Data Filtering and Extraction based on Study Area

# mexico_city <- municipalities_COVID_sf %>%
#   filter(CVE_ENT == "09")
# 
# mexico_state <- municipalities_COVID_sf %>%
#   filter(CVE_ENT == "15")
# 
# morelos_state <- municipalities_COVID_sf %>%
#   filter(CVE_ENT == "17")

mexico_city <- municipalities_COVID[(municipalities_COVID$CVE_ENT == '09'),]
mexico_state <- municipalities_COVID[(municipalities_COVID$CVE_ENT == '15'),]
morelos_state <- municipalities_COVID[(municipalities_COVID$CVE_ENT == '17'),]

mexi_combined <- municipalities_COVID[(municipalities_COVID$CVE_ENT == '09' | municipalities_COVID$CVE_ENT == '15' | municipalities_COVID$CVE_ENT == '17'),]
# 
# mexi_combined <- st_as_sf(mexi_combined)
# 
# mexi_combined <- mexi_combined %>%
#   filter("CVE_ENT" == "09") %>%
#     mutate("area" = "Mexico City")
# 
# mexi_combined <- mexi_combined %>%
#   filter("CVE_ENT" == "15") %>%
#   mutate("area" = "Mexico State")
# 
# mexi_combined <- mexi_combined %>%
#   filter("CVE_ENT" == "17") %>%
#   mutate("area" = "Morelos State")

3.2.1 Check for any NA

mexico_city <- na.omit(mexico_city)
mexico_state <- na.omit(mexico_state)
morelos_state <- na.omit(morelos_state)
mexi_combined <- na.omit(mexi_combined)

3.2.2 Plot the Study Area

CVE_ENT code which represent the area as follow: 09 - Mexico City 15 - Mexico State 17 - Morelos State

plot(mexico_city, main = "Study Area of Mexico City")

plot(mexico_state, main = "Study area of Mexico State")

plot(morelos_state, main = "Study Area of Morelos State")

plot(mexi_combined, main = "Combined of the Three Study Area")

qtm(mexi_combined, "CVE_ENT",  title = "Combined of the Three Study Area") +
  tm_legend(show = TRUE) + tm_text("NOMGEO", size=0.2) 

4 Exploratiory Data Analysis

Plot based on the data analyses on Cumulative, Active and Death rate of the Covid 19.

4.1 Filter Selection for Cumulative, Active and Death Cases

Add the sum of week 13 to 32 of new cases to achieve the cumulative cases Similarly we can perform the same method for active and death cases but using the respective column accordingly.

#c(1:6) --> dplyr::select("CVEGEO", "CVE_ENT", "CVE_MUN", "NOMGEO", "Pop2020", "wk13_32_total_count")

mexico_city_filtered <- mexico_city[,c(1:6, 19:38, 51:70, 83:102, 115:134)]
mexico_city_filtered_cumu <- mexico_city[,c(1:6,19:38)]
mexico_state_filtered <- mexico_state[,c(1:6, 19:38, 51:70, 83:102, 115:134)]
morelos_state_filtered <- morelos_state[,c(1:6, 19:38, 51:70, 83:102, 115:134)]
mexi_combined_filtered <- mexi_combined[,c(1:6, 19:38, 51:70, 83:102, 115:134)]

mexico_city_filtered@data$`wk13_32_cumul_count` <-  as.numeric(apply(mexico_city_filtered@data[,7:26], 1, sum))
#mexico_city_filtered@data$`wk13_32_cumulative` <-  mexico_city_filtered@data$`cumul32` - mexico_city_filtered@data$`cumul12`
mexico_city_filtered@data$`wk13_32_activ_count` <-  as.numeric(apply(mexico_city_filtered@data[,47:66], 1, sum))
mexico_city_filtered@data$`wk13_32_death_count` <-  as.numeric(apply(mexico_city_filtered@data[,67:86], 1, sum))

mexico_state_filtered@data$`wk13_32_cumul_count` <-  as.numeric(apply(mexico_state_filtered@data[,7:26], 1, sum))
#mexico_state_filtered@data$`wk13_32_cumulative` <-  mexico_state_filtered@data$`cumul32` - mexico_state_filtered@data$`cumul12`
mexico_state_filtered@data$`wk13_32_activ_count` <-  as.numeric(apply(mexico_state_filtered@data[,47:66], 1, sum))
mexico_state_filtered@data$`wk13_32_death_count` <-  as.numeric(apply(mexico_state_filtered@data[,67:86], 1, sum))

morelos_state_filtered@data$`wk13_32_cumul_count` <-  as.numeric(apply(morelos_state_filtered@data[,7:26], 1, sum))
#morelos_state_filtered@data$`wk13_32_cumulative` <-  morelos_state_filtered@data$`cumul32` - morelos_state_filtered@data$`cumul12`
morelos_state_filtered@data$`wk13_32_activ_count` <-  as.numeric(apply(morelos_state_filtered@data[,47:66], 1, sum))
morelos_state_filtered@data$`wk13_32_death_count` <-  as.numeric(apply(morelos_state_filtered@data[,67:86], 1, sum))

mexi_combined_filtered@data$`wk13_32_cumul_count` <-  as.numeric(apply(mexi_combined_filtered@data[,7:26], 1, sum))
#mexi_combined_filtered@data$`wk13_32_cumulative` <-  mexi_combined_filtered@data$`cumul32` - mexi_combined_filtered@data$`cumul12`
mexi_combined_filtered@data$`wk13_32_activ_count` <-  as.numeric(apply(mexi_combined_filtered@data[,47:66], 1, sum))
mexi_combined_filtered@data$`wk13_32_death_count` <-  as.numeric(apply(mexi_combined_filtered@data[,67:86], 1, sum))

4.2 Find out the rate of covid through Cumulative cases

For effective comparison, calculate the COVID-19 rate from e-week 13 till 32 based on cumulative cases Rate by per case per 10,000 population

mexico_city_filtered@data$`rate_of_covid` <- (mexico_city_filtered@data$`wk13_32_cumul_count` * 10000)/ mexico_city_filtered@data$`Pop2020`

mexico_state_filtered@data$`rate_of_covid` <- (mexico_state_filtered@data$`wk13_32_cumul_count` * 10000)/ mexico_state_filtered@data$`Pop2020`

morelos_state_filtered@data$`rate_of_covid` <- (morelos_state_filtered@data$`wk13_32_cumul_count` * 10000)/ morelos_state_filtered@data$`Pop2020`

mexi_combined_filtered@data$`rate_of_covid` <- (mexi_combined_filtered@data$`wk13_32_cumul_count` * 10000)/ mexi_combined_filtered@data$`Pop2020`

mexico_city_filtered_roc_sf <- st_as_sf(mexico_city_filtered) %>%
  select("CVE_ENT", "NOMGEO", "rate_of_covid")

mexico_state_filtered_roc_sf <- st_as_sf(mexico_state_filtered) %>%
  select("CVE_ENT", "NOMGEO", "rate_of_covid")

morelos_state_filtered_roc_sf <- st_as_sf(morelos_state_filtered) %>%
  select("CVE_ENT", "NOMGEO", "rate_of_covid")

mexi_combined_filtered_roc_sf <- st_as_sf(mexi_combined_filtered) %>%
  select("CVE_ENT", "NOMGEO", "rate_of_covid")

4.3 Find out the active rate of COVID 19

Calculating the COVID-19 rate from e-week 13 till 32 based on active cases Rate by per case per 10,000 population

mexico_city_filtered@data$`active_rate_of_covid` <- (mexico_city_filtered@data$`wk13_32_activ_count` * 10000)/ mexico_city_filtered@data$`Pop2020`

mexico_state_filtered@data$`active_rate_of_covid` <- (mexico_state_filtered@data$`wk13_32_activ_count` * 10000)/ mexico_state_filtered@data$`Pop2020`

morelos_state_filtered@data$`active_rate_of_covid` <- (morelos_state_filtered@data$`wk13_32_activ_count` * 10000)/ morelos_state_filtered@data$`Pop2020`

mexi_combined_filtered@data$`active_rate_of_covid` <- (mexi_combined_filtered@data$`wk13_32_activ_count` * 10000)/ mexi_combined_filtered@data$`Pop2020`

#print(mexico_city_filtered)
#print(mexico_state_filtered)
#print(morelos_state_filtered)
#print(mexi_combined_filtered)

4.4 Find out the death rate of COVID 19

Calculating the COVID-19 rate from e-week 13 till 32 based on death cases Rate by per case per 10,000 population

mexico_city_filtered@data$`death_rate_of_covid` <- (mexico_city_filtered@data$`wk13_32_death_count` * 10000)/ mexico_city_filtered@data$`Pop2020`

mexico_state_filtered@data$`death_rate_of_covid` <- (mexico_state_filtered@data$`wk13_32_death_count` * 10000)/ mexico_state_filtered@data$`Pop2020`

morelos_state_filtered@data$`death_rate_of_covid` <- (morelos_state_filtered@data$`wk13_32_death_count` * 10000)/ morelos_state_filtered@data$`Pop2020`

mexi_combined_filtered@data$`death_rate_of_covid` <- (mexi_combined_filtered@data$`wk13_32_death_count` * 10000)/ mexi_combined_filtered@data$`Pop2020`


mexico_city_filtered_sf <- st_as_sf(mexico_city_filtered)
mexico_state_filtered_sf <- st_as_sf(mexico_state_filtered)
morelos_state_filtered_sf <- st_as_sf(morelos_state_filtered)
mexi_combined_filtered_sf <- st_as_sf(mexi_combined_filtered)

#print(mexico_city_filtered_sf)
#print(mexico_state_filtered_sf)
#print(morelos_state_filtered_sf)
#print(mexi_combined_filtered_sf)

5 Exploratory Data Analysis

5.1 Barplots for Catergorical Variables

5.1.1 Cumulative COVID 19 Cases in each Municipalities for each City or States

#Size bottom left top right
par(mar=c(7, 5, 2, 1))
mexico_city_filtered_bp <- barplot(mexico_city_filtered$wk13_32_cumul_count, ylim = c(0, 16000), main = "Total No. of COVID 19 cases in each Municipalities of Mexico City", ylab="Cumulative Count")
axis(1, at=mexico_city_filtered_bp, labels=mexico_city_filtered$NOMGEO, las=3, cex.lab=4, cex.axis=0.7)

#Size bottom left top right
par(mar=c(7, 5, 2, 1))
mexico_state_filtered_bp <- barplot(mexico_state_filtered$wk13_32_cumul_count, ylim = c(0, 8000), main = "Total No. of COVID 19 cases in each Municipalities of Mexico State", ylab="Cumulative Count")
axis(1, at=mexico_state_filtered_bp, labels=mexico_state_filtered$NOMGEO, las=3, cex.lab=4, cex.axis=0.7)

#Size bottom left top right
par(mar=c(7, 5, 2, 1))
morelos_state_filtered_bp <- barplot(morelos_state_filtered$wk13_32_cumul_count, ylim = c(0, 1200), main = "Total No. of COVID 19 cases in each Municipalities of Morelos State", ylab="Cumulative Count")
axis(1, at=morelos_state_filtered_bp, labels=morelos_state_filtered$NOMGEO, las=3, cex.lab=4, cex.axis=0.7)

5.1.2 Rate of COVID 19 in each Municipalities for each City or States

#Size bottom left top right
par(mar=c(7, 5, 2, 1))
mexico_city_filtered_roc_bp <- barplot(mexico_city_filtered_roc_sf$rate_of_covid, ylim = c(0, 200), main = "Rate of COVID 19 in each Municipalities of Mexico City", ylab="Rate of COVID 19")
axis(1, at=mexico_city_filtered_roc_bp, labels=mexico_city_filtered_roc_sf$NOMGEO, las=3, cex.lab=4, cex.axis=0.7)

#Size bottom left top right
par(mar=c(7, 5, 2, 1))
mexico_state_filtered_roc_bp <- barplot(mexico_state_filtered_roc_sf$rate_of_covid, ylim = c(0, 200), main = "Rate of COVID 19 in each Municipalities of Mexico State", ylab="Rate of COVID 19")
axis(1, at=mexico_state_filtered_roc_bp, labels=mexico_state_filtered_roc_sf$NOMGEO, las=3, cex.lab=4, cex.axis=0.7)

#Size bottom left top right
par(mar=c(7, 5, 2, 1))
morelos_state_filtered_roc_bp <- barplot(morelos_state_filtered_roc_sf$rate_of_covid, ylim = c(0, 60), main = "Rate of COVID 19 in each Municipalities of Morelos State", ylab="Rate of COVID 19")
axis(1, at=morelos_state_filtered_roc_bp, labels=morelos_state_filtered_roc_sf$NOMGEO, las=3, cex.lab=4, cex.axis=0.7)

## Histogram to compare Frequency for Municipalities of population and COVID 19 Rate

hist(mexi_combined$Pop2020, main = "Municipality Population Distribution")

hist(mexi_combined_filtered$`rate_of_covid`, main = "Municipality COVID Rate Distribution (Per 10,000)")

5.2 Choloropleth Mapping

COVID 19 Rate of Mexico City (Per 10,000)

tm_shape(mexico_city_filtered) +
  tm_fill("rate_of_covid",
          style = "quantile", 
          thres.poly = 0,
          n = 5,
          palette = "Oranges",
          title = "COVID 19 Rate of Mexico City (Per 10,000)") +
 tm_layout(legend.show = TRUE, title.position = c("center", "center"),
            title.size = 20) +
  tm_borders(alpha = 0.5) + tm_text("NOMGEO", size=0.2) 

COVID 19 Rate of Mexico State (Per 10,000)

tm_shape(mexico_state_filtered) +
  tm_fill("rate_of_covid",
          style = "quantile", 
          thres.poly = 0,
          n = 5,
          palette = "Oranges",
          title = "COVID 19 Rate of Mexico State (Per 10,000)") +
 tm_layout(legend.show = TRUE, title.position = c("center", "center"),
            title.size = 20) +
  tm_borders(alpha = 0.5) + tm_text("NOMGEO", size=0.2) 

COVID 19 Rate of Morelos State (Per 10,000)

tm_shape(morelos_state_filtered) +
  tm_fill("rate_of_covid",
          style = "quantile", 
          thres.poly = 0,
          n = 5,
          palette = "Oranges",
          title = "COVID 19 Rate of Morelos State (Per 10,000)") +
 tm_layout(legend.show = TRUE, title.position = c("center", "center"),
            title.size = 20) +
  tm_borders(alpha = 0.5) + tm_text("NOMGEO", size=0.2) 

COVID 19 Rate of Three Selected Study Area (Per 10,000)

tm_shape(mexi_combined_filtered) +
  tm_fill("rate_of_covid",
          style = "quantile", 
          thres.poly = 0,
          n = 5,
          palette = "Oranges",
          title = "COVID 19 Rate of Three Selected Study Area (Per 10,000)") +
 tm_layout(legend.show = TRUE, title.position = c("center", "center"),
            title.size = 20) +
  tm_borders(alpha = 0.5) + tm_text("NOMGEO", size=0.2) 

Based on the above result, Mexico City have the most highest concentration of COVID 19 Rate.

5.3 Box Map Method

Creating Boxbreaks function

boxbreaks <- function(v, mult=1.5) { 
  qv <- unname(quantile(v))
  iqr <- qv[4] - qv[2]
  upfence <- qv[4] + mult * iqr 
  lofence <- qv[2] - mult * iqr
  bb <- vector(mode="numeric",length=7)

  if (lofence < qv[1]) { # no lower outliers
    bb[1] <- lofence
    bb[2] <- floor(qv[1]) 
    } 
  else {
    bb[2] <- lofence
    bb[1] <- qv[1] 
    }
  if (upfence > qv[5]) { # no upper outliers
    bb[7] <- upfence
    bb[6] <- ceiling(qv[5])
    } 
  else {
    bb[6] <- upfence
    bb[7] <- qv[5]
    }
  bb[3:5] <- qv[2:4]
  return(bb)
}

Creating the get.var function

get.var <- function(vname,df) { 
  v <- df[vname] %>% st_set_geometry(NULL)
  v <- unname(v[,1])
  return(v)
}

Box Map Function

boxmap <- function(vnam, df, legtitle=NA, mtitle="Box Map", mult=1.5){
  var <- get.var(vnam,df)
  bb <- boxbreaks(var)
  tm_shape(df) + 
  tm_fill(vnam, title=legtitle, breaks=bb, palette="-RdBu", labels = c("lower outlier", "< 25%", "25% - 50%", "50% - 75%","> 75%", "upper outlier")) +
  tm_borders() + 
  tm_layout(main.title = "Box Map of COVID-19 Rate in our Study Area", main.title.position = "center", main.title.size = 1.3, legend.position = c("right", "top"), legend.outside = TRUE , frame = FALSE)
}  + tm_text("NOMGEO", size=0.15) 

Mapping out the result

boxmap("rate_of_covid", mexi_combined_filtered_sf)

The result above allows us to better identify the Outlier in the Mexico city.

6 Modelling of the Spatial Neighbours

Using poly2nb() of spdep package to compute contiguity weight matrices for the study area. This function builds a neighbours list based on regions with contiguous boundaries.

6.1 Creating (QUEEN) contiguity based neighbours

wm_q_mexi_combined_filtered <- poly2nb(mexi_combined_filtered, queen=TRUE)
summary(wm_q_mexi_combined_filtered)
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 962 
## Percentage nonzero weights: 3.10563 
## Average number of links: 5.465909 
## Link number distribution:
## 
##  1  2  3  4  5  6  7  8  9 10 11 14 
##  3  3 18 41 31 31 22 13  9  3  1  1 
## 3 least connected regions:
## 721 739 765 with 1 link
## 1 most connected region:
## 761 with 14 links

The summary report above shows that there are 176 area units in Mexico City, Mexico State and Morelos State. The most connect area unit has 14 neighbours. There are 1 area unit with only 11 neighbours.

List complete weight matrix

str(wm_q_mexi_combined_filtered)
## List of 176
##  $ : int [1:5] 4 14 15 73 120
##  $ : int [1:5] 6 9 11 12 13
##  $ : int [1:4] 9 15 53 78
##  $ : int [1:8] 1 14 16 36 49 74 120 125
##  $ : int [1:5] 6 13 14 16 74
##  $ : int [1:8] 2 5 10 12 13 74 86 138
##  $ : int [1:3] 9 11 78
##  $ : int [1:8] 10 11 12 41 66 150 161 164
##  $ : int [1:7] 2 3 7 11 13 15 78
##  $ : int [1:5] 6 8 12 41 138
##  $ : int [1:9] 2 7 8 9 12 59 78 117 150
##  $ : int [1:5] 2 6 8 10 11
##  $ : int [1:6] 2 5 6 9 14 15
##  $ : int [1:6] 1 4 5 13 15 16
##  $ : int [1:7] 1 3 9 13 14 53 73
##  $ : int [1:4] 4 5 14 74
##  $ : int [1:4] 19 30 101 118
##  $ : int [1:7] 27 44 49 97 108 109 116
##  $ : int [1:4] 17 61 87 118
##  $ : int [1:4] 37 96 113 133
##  $ : int [1:7] 23 58 90 103 122 130 134
##  $ : int [1:4] 28 89 114 117
##  $ : int [1:7] 21 48 102 126 127 130 134
##  $ : int [1:3] 96 98 121
##  $ : int [1:4] 31 33 84 119
##  $ : int [1:2] 52 112
##  $ : int [1:6] 18 44 46 49 115 116
##  $ : int [1:3] 22 89 117
##  $ : int [1:6] 54 62 73 76 120 137
##  $ : int [1:5] 17 64 72 101 118
##  $ : int [1:5] 25 50 84 163 171
##  $ : int [1:4] 77 81 91 100
##  $ : int [1:6] 25 66 84 105 110 119
##  $ : int [1:8] 43 70 71 88 89 106 117 122
##  $ : int [1:3] 67 78 117
##  $ : int [1:4] 4 49 124 125
##  $ : int [1:8] 20 56 102 106 113 129 133 134
##  $ : int [1:3] 41 99 119
##  $ : int [1:4] 51 107 111 136
##  $ : int [1:8] 69 75 107 111 124 125 136 137
##  $ : int [1:9] 8 10 38 55 66 99 105 119 138
##  $ : int [1:4] 61 72 118 128
##  $ : int [1:4] 34 70 71 117
##  $ : int [1:7] 18 27 46 85 109 115 116
##  $ : int [1:4] 47 55 86 115
##  $ : int [1:3] 27 44 115
##  $ : int [1:4] 45 74 86 115
##  $ : int [1:4] 23 57 126 127
##  $ : int [1:11] 4 18 27 36 60 74 97 115 120 125 ...
##  $ : int [1:4] 31 157 163 171
##  $ : int [1:4] 39 111 112 136
##  $ : int [1:3] 26 112 136
##  $ : int [1:5] 3 15 67 73 78
##  $ : int [1:5] 29 62 76 83 103
##  $ : int [1:6] 41 45 86 115 119 138
##  $ : int [1:5] 37 123 129 133 135
##  $ : int [1:3] 48 94 126
##  $ : int [1:5] 21 63 64 90 103
##  $ : int [1:4] 11 79 117 150
##  $ : int [1:6] 49 75 97 125 136 141
##  $ : int [1:5] 19 42 87 95 118
##  $ : int [1:5] 29 54 73 83 131
##  $ : int [1:6] 58 64 72 76 103 128
##  $ : int [1:7] 30 58 63 72 80 90 101
##  $ : int [1:6] 68 79 104 106 114 117
##  $ : int [1:9] 8 33 41 84 99 105 110 164 168
##  $ : int [1:10] 35 53 70 73 78 83 92 117 122 131
##  $ : int [1:6] 65 79 104 135 146 156
##  $ : int [1:3] 40 75 124
##  $ : int [1:7] 34 43 67 71 92 117 122
##  $ : int [1:3] 34 43 70
##  $ : int [1:6] 30 42 63 64 118 128
##  $ : int [1:8] 1 15 29 53 62 67 120 131
##  $ : int [1:8] 4 5 6 16 47 49 86 115
##  $ : int [1:8] 40 60 69 97 124 125 136 141
##  $ : int [1:7] 29 54 63 103 111 128 137
##  $ : int 32
##  $ : int [1:8] 3 7 9 11 35 53 67 117
##  $ : int [1:8] 59 65 68 104 117 148 150 156
##  $ : int [1:4] 64 90 101 140
##  $ : int [1:4] 32 91 108 109
##  $ : int [1:5] 94 102 126 132 139
##  $ : int [1:6] 54 62 67 103 122 131
##  $ : int [1:7] 25 31 33 66 110 143 171
##  $ : int [1:3] 44 109 115
##  $ : int [1:6] 6 45 47 55 74 138
##  $ : int [1:2] 19 61
##  $ : int [1:4] 34 89 106 114
##  $ : int [1:6] 22 28 34 88 114 117
##  $ : int [1:6] 21 58 64 80 130 140
##  $ : int [1:5] 32 81 100 108 109
##  $ : int [1:3] 67 70 122
##  $ : int [1:3] 98 102 113
##  $ : int [1:3] 57 82 126
##  $ : int 61
##  $ : int [1:5] 20 24 98 113 133
##  $ : int [1:8] 18 49 60 75 100 108 136 141
##  $ : int [1:7] 24 93 96 102 113 132 139
##  $ : int [1:5] 38 41 66 105 119
##   [list output truncated]
##  - attr(*, "class")= chr "nb"
##  - attr(*, "region.id")= chr [1:176] "270" "271" "272" "273" ...
##  - attr(*, "call")= language poly2nb(pl = mexi_combined_filtered, queen = TRUE)
##  - attr(*, "type")= chr "queen"
##  - attr(*, "sym")= logi TRUE

6.1.1 Creating (ROOK) contiguity based neighbours

wm_r_mexi_combined_filtered <- poly2nb(mexi_combined_filtered, queen=FALSE)
summary(wm_r_mexi_combined_filtered)
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 946 
## Percentage nonzero weights: 3.053977 
## Average number of links: 5.375 
## Link number distribution:
## 
##  1  2  3  4  5  6  7  8  9 10 14 
##  3  3 18 47 27 31 22 14  6  4  1 
## 3 least connected regions:
## 721 739 765 with 1 link
## 1 most connected region:
## 761 with 14 links

The summary report above shows that there are 176 area units in Mexico City, Mexico State and Morelos State. The most connect area unit has 14 neighbours. There are 4 area units with 10 neighbours.

List complete weight matrix

str(wm_r_mexi_combined_filtered)
## List of 176
##  $ : int [1:5] 4 14 15 73 120
##  $ : int [1:5] 6 9 11 12 13
##  $ : int [1:4] 9 15 53 78
##  $ : int [1:8] 1 14 16 36 49 74 120 125
##  $ : int [1:5] 6 13 14 16 74
##  $ : int [1:8] 2 5 10 12 13 74 86 138
##  $ : int [1:3] 9 11 78
##  $ : int [1:8] 10 11 12 41 66 150 161 164
##  $ : int [1:7] 2 3 7 11 13 15 78
##  $ : int [1:5] 6 8 12 41 138
##  $ : int [1:9] 2 7 8 9 12 59 78 117 150
##  $ : int [1:5] 2 6 8 10 11
##  $ : int [1:6] 2 5 6 9 14 15
##  $ : int [1:6] 1 4 5 13 15 16
##  $ : int [1:7] 1 3 9 13 14 53 73
##  $ : int [1:4] 4 5 14 74
##  $ : int [1:4] 19 30 101 118
##  $ : int [1:7] 27 44 49 97 108 109 116
##  $ : int [1:4] 17 61 87 118
##  $ : int [1:4] 37 96 113 133
##  $ : int [1:7] 23 58 90 103 122 130 134
##  $ : int [1:4] 28 89 114 117
##  $ : int [1:7] 21 48 102 126 127 130 134
##  $ : int [1:3] 96 98 121
##  $ : int [1:4] 31 33 84 119
##  $ : int [1:2] 52 112
##  $ : int [1:6] 18 44 46 49 115 116
##  $ : int [1:3] 22 89 117
##  $ : int [1:6] 54 62 73 76 120 137
##  $ : int [1:5] 17 64 72 101 118
##  $ : int [1:5] 25 50 84 163 171
##  $ : int [1:4] 77 81 91 100
##  $ : int [1:5] 25 66 84 105 119
##  $ : int [1:8] 43 70 71 88 89 106 117 122
##  $ : int [1:3] 67 78 117
##  $ : int [1:4] 4 49 124 125
##  $ : int [1:8] 20 56 102 106 113 129 133 134
##  $ : int [1:3] 41 99 119
##  $ : int [1:4] 51 107 111 136
##  $ : int [1:8] 69 75 107 111 124 125 136 137
##  $ : int [1:8] 8 10 38 55 66 99 119 138
##  $ : int [1:4] 61 72 118 128
##  $ : int [1:4] 34 70 71 117
##  $ : int [1:7] 18 27 46 85 109 115 116
##  $ : int [1:4] 47 55 86 115
##  $ : int [1:3] 27 44 115
##  $ : int [1:4] 45 74 86 115
##  $ : int [1:4] 23 57 126 127
##  $ : int [1:10] 4 18 27 36 60 74 97 115 120 125
##  $ : int [1:4] 31 157 163 171
##  $ : int [1:4] 39 111 112 136
##  $ : int [1:3] 26 112 136
##  $ : int [1:5] 3 15 67 73 78
##  $ : int [1:5] 29 62 76 83 103
##  $ : int [1:6] 41 45 86 115 119 138
##  $ : int [1:5] 37 123 129 133 135
##  $ : int [1:3] 48 94 126
##  $ : int [1:5] 21 63 64 90 103
##  $ : int [1:4] 11 79 117 150
##  $ : int [1:5] 49 75 125 136 141
##  $ : int [1:5] 19 42 87 95 118
##  $ : int [1:5] 29 54 73 83 131
##  $ : int [1:6] 58 64 72 76 103 128
##  $ : int [1:7] 30 58 63 72 80 90 101
##  $ : int [1:6] 68 79 104 106 114 117
##  $ : int [1:7] 8 33 41 105 110 164 168
##  $ : int [1:10] 35 53 70 73 78 83 92 117 122 131
##  $ : int [1:6] 65 79 104 135 146 156
##  $ : int [1:3] 40 75 124
##  $ : int [1:7] 34 43 67 71 92 117 122
##  $ : int [1:3] 34 43 70
##  $ : int [1:6] 30 42 63 64 118 128
##  $ : int [1:8] 1 15 29 53 62 67 120 131
##  $ : int [1:8] 4 5 6 16 47 49 86 115
##  $ : int [1:8] 40 60 69 97 124 125 136 141
##  $ : int [1:7] 29 54 63 103 111 128 137
##  $ : int 32
##  $ : int [1:8] 3 7 9 11 35 53 67 117
##  $ : int [1:8] 59 65 68 104 117 148 150 156
##  $ : int [1:4] 64 90 101 140
##  $ : int [1:4] 32 91 108 109
##  $ : int [1:4] 94 126 132 139
##  $ : int [1:6] 54 62 67 103 122 131
##  $ : int [1:6] 25 31 33 110 143 171
##  $ : int [1:3] 44 109 115
##  $ : int [1:6] 6 45 47 55 74 138
##  $ : int [1:2] 19 61
##  $ : int [1:4] 34 89 106 114
##  $ : int [1:6] 22 28 34 88 114 117
##  $ : int [1:6] 21 58 64 80 130 140
##  $ : int [1:5] 32 81 100 108 109
##  $ : int [1:3] 67 70 122
##  $ : int [1:3] 98 102 113
##  $ : int [1:3] 57 82 126
##  $ : int 61
##  $ : int [1:5] 20 24 98 113 133
##  $ : int [1:7] 18 49 75 100 108 136 141
##  $ : int [1:7] 24 93 96 102 113 132 139
##  $ : int [1:4] 38 41 105 119
##   [list output truncated]
##  - attr(*, "class")= chr "nb"
##  - attr(*, "region.id")= chr [1:176] "270" "271" "272" "273" ...
##  - attr(*, "call")= language poly2nb(pl = mexi_combined_filtered, queen = FALSE)
##  - attr(*, "type")= chr "rook"
##  - attr(*, "sym")= logi TRUE

6.1.2 Plotting Queen contiguity against ROOK contiguity based neighbours maps

par(mfrow=c(1,2))
plot(mexi_combined_filtered, border="lightgrey")
plot(wm_q_mexi_combined_filtered, coordinates(mexi_combined_filtered), pch = 19, cex = 0.6, add = TRUE, col = "red", ain="Queen Contiguity")
plot(mexi_combined_filtered, border="lightgrey")
plot(wm_r_mexi_combined_filtered, coordinates(mexi_combined_filtered), pch = 19, cex = 0.6, add = TRUE, col = "red", ain="ROOK Contiguity")

## Computing distance based neighbours

To derive distance-based weight matrices by using dnearneigh() of spdep package.

The function identifies neighbours of region points by Euclidean distance with a distance band with lower d1= and upper d2= bounds controlled by the bounds= argument. If unprojected coordinates are used and either specified in the coordinates object x or with x as a two column matrix and longlat=TRUE, great circle distances in km will be calculated assuming the WGS84 reference ellipsoid.

6.1.3 Determine the cut-off distance

  • Return a matrix with the indices of points belonging to the set of the k nearest neighbours of each other by using knearneigh() of spdep.
  • Convert the knn object returned by knearneigh() into a neighbours list of class nb with a list of integer vectors containing neighbour region number ids by using knn2nb().
  • Return the length of neighbour relationship edges by using nbdists() of spdep. The function returns in the units of the coordinates if the coordinates are projected, in km otherwise.
  • Remove the list structure of the returned object by using unlist().
coords <- coordinates(mexi_combined_filtered)
k1 <- knn2nb(knearneigh(coords))
k1dists <- unlist(nbdists(k1, coords, longlat = TRUE))
summary(k1dists)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1267    6953   10703   10470   13830   19827

The summary report shows that the largest first nearest neighbour distance is 19827 km, so using this as the upper threshold gives certainty that all units will have at least one neighbour.

6.1.4 Computing fixed distance weight matrix

Now, we will compute the distance weight matrix by using dnearneigh() as shown in the code chunk below.

wm_d19827 <- dnearneigh(coords, 0, 19827, longlat = TRUE)
wm_d19827
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 30784 
## Percentage nonzero weights: 99.38017 
## Average number of links: 174.9091

Next, we will use str() to display the content of wm_d62 weight matrix.

str(wm_d19827)
## List of 176
##  $ : int [1:175] 2 3 4 5 6 7 8 9 10 11 ...
##  $ : int [1:175] 1 3 4 5 6 7 8 9 10 11 ...
##  $ : int [1:173] 1 2 4 5 6 7 8 9 10 11 ...
##  $ : int [1:175] 1 2 3 5 6 7 8 9 10 11 ...
##  $ : int [1:174] 1 2 3 4 6 7 8 10 11 12 ...
##  $ : int [1:175] 1 2 3 4 5 7 8 9 10 11 ...
##  $ : int [1:175] 1 2 3 4 5 6 8 9 10 11 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 9 10 11 ...
##  $ : int [1:174] 1 2 3 4 6 7 8 10 11 12 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 11 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:174] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:174] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:174] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:174] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:174] 1 2 3 4 5 6 7 8 9 10 ...
##  $ : int [1:175] 1 2 3 4 5 6 7 8 9 10 ...
##   [list output truncated]
##  - attr(*, "class")= chr "nb"
##  - attr(*, "nbtype")= chr "distance"
##  - attr(*, "distances")= num [1:2] 0 19827
##  - attr(*, "region.id")= chr [1:176] "1" "2" "3" "4" ...
##  - attr(*, "call")= language dnearneigh(x = coords, d1 = 0, d2 = 19827, longlat = TRUE)
##  - attr(*, "dnn")= num [1:2] 0 19827
##  - attr(*, "bounds")= chr [1:2] "GT" "LE"
##  - attr(*, "sym")= logi TRUE

Using table() and card() of spdep to display the structure of the weight matrix

table(mexi_combined_filtered$CVE_ENT, card(wm_d19827))
##     
##      173 174 175
##   09   1   2  13
##   15   0   8 117
##   17   1   2  32
n_comp <- n.comp.nb(wm_d19827)
n_comp$nc
## [1] 1
table(n_comp$comp.id)
## 
##   1 
## 176

6.1.5 Plotting fixed distance weight matrix

plot(mexi_combined_filtered, border="lightgrey")
plot(wm_d19827, coords, add=TRUE)
plot(k1, coords, add=TRUE, col="red", length=0.08)

The red lines show the links of 1st nearest neighbours and the black lines show the links of neighbours within the cut-off distance of 19827km.

par(mfrow=c(1,2))
plot(mexi_combined_filtered, border="lightgrey")
plot(k1, coords, add=TRUE, col="red", length=0.08, main="1st nearest neighbours")
plot(mexi_combined_filtered, border="lightgrey")
plot(wm_d19827, coords, add=TRUE, pch = 19, cex = 0.6, main="Distance link")

### Computing adaptive distance weight matrix

One of the characteristics of fixed distance weight matrix is that more densely settled areas (usually the urban areas) tend to have more neighbours and the less densely settled areas (usually the rural counties) tend to have lesser neighbours. Having many neighbours smoothes the neighbour relationship across more neighbours.

It is possible to control the numbers of neighbours directly using k-nearest neighbours, either accepting asymmetric neighbours or imposing symmetry as shown in the code chunk below.

knn6 <- knn2nb(knearneigh(coords, k=6))
knn6
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 1056 
## Percentage nonzero weights: 3.409091 
## Average number of links: 6 
## Non-symmetric neighbours list
str(knn6)
## List of 176
##  $ : int [1:6] 4 14 15 16 73 120
##  $ : int [1:6] 5 6 9 12 13 14
##  $ : int [1:6] 7 9 15 53 73 78
##  $ : int [1:6] 1 5 14 16 49 120
##  $ : int [1:6] 2 6 13 14 16 74
##  $ : int [1:6] 2 5 10 13 16 74
##  $ : int [1:6] 2 3 9 11 53 78
##  $ : int [1:6] 10 11 12 66 161 164
##  $ : int [1:6] 2 3 7 13 15 53
##  $ : int [1:6] 6 8 12 86 99 138
##  $ : int [1:6] 2 7 9 12 59 150
##  $ : int [1:6] 2 6 8 10 11 138
##  $ : int [1:6] 2 5 9 14 15 16
##  $ : int [1:6] 1 4 5 13 15 16
##  $ : int [1:6] 1 9 13 14 16 73
##  $ : int [1:6] 4 5 6 13 14 74
##  $ : int [1:6] 19 30 64 87 101 118
##  $ : int [1:6] 44 46 85 97 108 116
##  $ : int [1:6] 17 30 61 87 101 118
##  $ : int [1:6] 37 56 93 96 113 133
##  $ : int [1:6] 23 58 103 122 130 134
##  $ : int [1:6] 28 35 43 89 114 117
##  $ : int [1:6] 48 102 126 127 130 134
##  $ : int [1:6] 93 96 98 121 132 139
##  $ : int [1:6] 31 33 50 84 105 119
##  $ : int [1:6] 39 51 52 107 112 136
##  $ : int [1:6] 18 44 46 49 85 116
##  $ : int [1:6] 22 35 43 89 114 117
##  $ : int [1:6] 1 62 73 120 125 137
##  $ : int [1:6] 17 64 72 90 101 118
##  $ : int [1:6] 25 33 50 84 110 163
##  $ : int [1:6] 77 81 91 100 108 109
##  $ : int [1:6] 25 38 66 84 99 105
##  $ : int [1:6] 43 70 71 88 89 106
##  $ : int [1:6] 22 28 43 78 114 117
##  $ : int [1:6] 40 49 69 124 125 141
##  $ : int [1:6] 20 56 104 106 113 129
##  $ : int [1:6] 33 41 99 105 119 138
##  $ : int [1:6] 40 51 60 69 107 111
##  $ : int [1:6] 60 69 75 107 124 137
##  $ : int [1:6] 38 55 99 105 119 138
##  $ : int [1:6] 61 72 76 95 118 128
##  $ : int [1:6] 28 70 71 88 89 92
##  $ : int [1:6] 18 27 46 85 109 116
##  $ : int [1:6] 47 55 74 86 115 138
##  $ : int [1:6] 18 27 44 85 115 116
##  $ : int [1:6] 6 27 45 74 86 115
##  $ : int [1:6] 23 57 94 126 127 130
##  $ : int [1:6] 4 27 36 116 125 141
##  $ : int [1:6] 31 84 110 143 157 163
##  $ : int [1:6] 39 40 69 107 111 112
##  $ : int [1:6] 26 60 75 100 112 136
##  $ : int [1:6] 3 7 9 62 73 131
##  $ : int [1:6] 29 62 73 76 83 131
##  $ : int [1:6] 38 41 45 86 115 119
##  $ : int [1:6] 20 37 104 123 129 135
##  $ : int [1:6] 48 82 94 126 127 132
##  $ : int [1:6] 21 30 63 64 90 103
##  $ : int [1:6] 7 22 28 35 78 117
##  $ : int [1:6] 40 69 75 107 124 136
##  $ : int [1:6] 17 19 42 72 95 118
##  $ : int [1:6] 29 53 54 73 76 131
##  $ : int [1:6] 54 58 72 83 103 128
##  $ : int [1:6] 30 58 63 72 90 118
##  $ : int [1:6] 22 28 88 89 104 114
##  $ : int [1:6] 33 99 105 110 164 168
##  $ : int [1:6] 35 53 78 83 92 131
##  $ : int [1:6] 79 104 123 135 146 156
##  $ : int [1:6] 40 60 75 107 124 141
##  $ : int [1:6] 34 43 71 89 92 122
##  $ : int [1:6] 34 43 70 88 89 92
##  $ : int [1:6] 30 42 63 64 118 128
##  $ : int [1:6] 1 15 29 53 62 120
##  $ : int [1:6] 4 5 6 16 47 86
##  $ : int [1:6] 40 60 69 124 136 141
##  $ : int [1:6] 29 54 62 111 128 137
##  $ : int [1:6] 32 81 91 100 108 109
##  $ : int [1:6] 3 28 35 59 67 117
##  $ : int [1:6] 65 68 114 117 148 150
##  $ : int [1:6] 17 30 64 90 101 140
##  $ : int [1:6] 32 77 85 91 108 109
##  $ : int [1:6] 48 57 94 126 132 139
##  $ : int [1:6] 54 62 63 67 103 131
##  $ : int [1:6] 31 33 50 66 110 168
##  $ : int [1:6] 18 44 46 109 115 116
##  $ : int [1:6] 6 10 45 47 74 138
##  $ : int [1:6] 17 19 61 95 101 118
##  $ : int [1:6] 22 28 43 71 89 114
##  $ : int [1:6] 22 28 43 71 88 114
##  $ : int [1:6] 30 58 64 80 130 140
##  $ : int [1:6] 18 32 77 81 100 108
##  $ : int [1:6] 35 43 67 70 71 89
##  $ : int [1:6] 20 98 102 113 126 132
##  $ : int [1:6] 48 57 82 126 127 132
##  $ : int [1:6] 19 42 61 72 118 128
##  $ : int [1:6] 20 24 93 98 113 133
##  $ : int [1:6] 18 60 75 108 116 141
##  $ : int [1:6] 24 82 93 96 132 139
##  $ : int [1:6] 33 38 41 66 105 138
##   [list output truncated]
##  - attr(*, "region.id")= chr [1:176] "1" "2" "3" "4" ...
##  - attr(*, "call")= language knearneigh(x = coords, k = 6)
##  - attr(*, "sym")= logi FALSE
##  - attr(*, "type")= chr "knn"
##  - attr(*, "knn-k")= num 6
##  - attr(*, "class")= chr "nb"

6.1.6 Plotting distance based neighbours

plot(mexi_combined_filtered, border="lightgrey")
plot(knn6, coords, pch = 19, cex = 0.6, add = TRUE, col = "red")

6.2 Weights based on IDW

Derive a spatial weight matrix based on Inversed Distance method. First, we will compute the distances between areas by using nbdists() of spdep.

dist <- nbdists(wm_q_mexi_combined_filtered, coords, longlat = TRUE)
ids <- lapply(dist, function(x) 1/(x))
ids
## [[1]]
## [1] 1.625282e-04 7.286620e-05 7.107408e-05 3.382463e-04 7.435495e-05
## 
## [[2]]
## [1] 5.811475e-05 6.758899e-05 5.642515e-05 7.017576e-05 8.656476e-05
## 
## [[3]]
## [1] 7.892491e-04 6.549298e-05 1.741204e-04 6.317576e-05
## 
## [[4]]
## [1] 1.625282e-04 8.303279e-05 6.389969e-05 7.233752e-05 1.165162e-04
## [6] 2.448648e-04 7.937698e-05 9.351276e-05
## 
## [[5]]
## [1] 5.822100e-05 6.080187e-05 1.057847e-04 8.574530e-05 9.734377e-05
## 
## [[6]]
## [1] 5.811475e-05 5.822100e-05 1.084275e-04 2.913919e-04 1.575909e-04
## [6] 1.447864e-04 7.795105e-05 7.082767e-05
## 
## [[7]]
## [1] 9.168834e-05 1.253962e-04 7.378774e-05
## 
## [[8]]
## [1] 8.051543e-05 9.717618e-05 1.392697e-04 8.476316e-05 1.178655e-04
## [6] 6.381720e-05 1.259769e-04 1.883358e-04
## 
## [[9]]
## [1] 6.758899e-05 7.892491e-04 9.168834e-05 3.408606e-04 2.768554e-04
## [6] 6.487048e-05 6.829824e-05
## 
## [[10]]
## [1] 1.084275e-04 8.051543e-05 9.885437e-05 2.440929e-04 1.367894e-04
## 
## [[11]]
## [1] 5.642515e-05 1.253962e-04 9.717618e-05 3.408606e-04 2.877919e-04
## [6] 2.639822e-04 5.920533e-05 6.716088e-05 7.704263e-05
## 
## [[12]]
## [1] 7.017576e-05 2.913919e-04 1.392697e-04 9.885437e-05 2.877919e-04
## 
## [[13]]
## [1] 8.656476e-05 6.080187e-05 1.575909e-04 2.768554e-04 8.288005e-05
## [6] 8.457247e-05
## 
## [[14]]
## [1] 7.286620e-05 8.303279e-05 1.057847e-04 8.288005e-05 8.401945e-05
## [6] 8.893802e-05
## 
## [[15]]
## [1] 7.107408e-05 6.549298e-05 6.487048e-05 8.457247e-05 8.401945e-05
## [6] 5.500956e-05 8.881606e-05
## 
## [[16]]
## [1] 6.389969e-05 8.574530e-05 8.893802e-05 5.586925e-05
## 
## [[17]]
## [1] 7.073320e-04 1.435673e-04 5.669883e-05 1.063840e-04
## 
## [[18]]
## [1] 1.356883e-04 7.528173e-05 1.276375e-04 9.455630e-05 6.162470e-05
## [6] 7.825048e-05 7.484700e-05
## 
## [[19]]
## [1] 7.073320e-04 5.074741e-05 1.147090e-04 1.232539e-04
## 
## [[20]]
## [1] 6.207069e-05 1.989998e-04 2.572725e-04 1.955575e-04
## 
## [[21]]
## [1] 6.151062e-05 1.274074e-04 1.794158e-04 9.416121e-05 1.114368e-04
## [6] 7.500310e-05 1.043238e-04
## 
## [[22]]
## [1] 7.501270e-05 7.591507e-05 2.928491e-04 2.055733e-04
## 
## [[23]]
## [1] 6.151062e-05 1.547248e-04 7.573288e-05 6.228594e-05 1.206285e-04
## [6] 1.150039e-04 8.822855e-05
## 
## [[24]]
## [1] 9.839306e-05 6.215702e-05 7.026734e-05
## 
## [[25]]
## [1] 1.772273e-04 1.597936e-04 7.456176e-05 7.246207e-05
## 
## [[26]]
## [1] 1.728270e-04 5.240886e-05
## 
## [[27]]
## [1] 1.356883e-04 5.745815e-05 9.279165e-05 8.762809e-05 1.489506e-04
## [6] 6.348349e-05
## 
## [[28]]
## [1] 7.501270e-05 8.074777e-05 1.046676e-04
## 
## [[29]]
## [1] 1.092489e-04 1.355575e-04 1.246360e-04 6.797692e-05 1.155047e-04
## [6] 6.979176e-05
## 
## [[30]]
## [1] 1.435673e-04 9.289630e-05 1.570112e-04 7.459317e-05 6.827966e-05
## 
## [[31]]
## [1] 1.772273e-04 1.016815e-04 5.614833e-05 3.042400e-04 9.680763e-05
## 
## [[32]]
## [1] 1.918575e-04 6.965538e-05 1.148537e-04 1.461145e-04
## 
## [[33]]
## [1] 1.597936e-04 1.126621e-04 1.142536e-04 2.320960e-04 3.757265e-03
## [6] 5.013834e-05
## 
## [[34]]
## [1] 5.359027e-04 5.982625e-05 6.981223e-05 7.253707e-05 6.435784e-05
## [6] 8.143763e-05 9.297839e-05 1.617230e-04
## 
## [[35]]
## [1] 6.333917e-05 8.473366e-05 1.036850e-04
## 
## [[36]]
## [1] 7.233752e-05 6.065400e-05 7.768842e-05 1.265698e-04
## 
## [[37]]
## [1] 6.207069e-05 6.851713e-05 1.162834e-04 8.780240e-05 6.920613e-05
## [6] 1.108079e-04 6.175376e-05 6.794225e-05
## 
## [[38]]
## [1] 7.531720e-05 5.348017e-05 5.986237e-05
## 
## [[39]]
## [1] 1.161859e-04 1.205144e-04 8.003154e-05 1.749139e-04
## 
## [[40]]
## [1] 6.310945e-05 1.502051e-04 5.043640e-05 6.274283e-05 6.583885e-05
## [6] 6.063476e-05 5.785814e-05 1.370508e-04
## 
## [[41]]
## [1] 8.476316e-05 2.440929e-04 7.531720e-05 1.414978e-04 5.261504e-05
## [6] 1.544491e-04 1.341811e-04 1.053480e-04 9.422543e-05
## 
## [[42]]
## [1] 2.958042e-04 7.717293e-05 1.130843e-04 1.658591e-04
## 
## [[43]]
## [1] 5.359027e-04 6.367085e-05 6.277917e-05 9.400085e-05
## 
## [[44]]
## [1] 7.528173e-05 5.745815e-05 9.571071e-05 7.616479e-05 1.400602e-04
## [6] 8.314604e-05 6.054047e-04
## 
## [[45]]
## [1] 1.473610e-04 1.292529e-04 1.071633e-04 6.943647e-05
## 
## [[46]]
## [1] 9.279165e-05 9.571071e-05 5.725829e-05
## 
## [[47]]
## [1] 1.473610e-04 5.328770e-05 7.801567e-05 7.605560e-05
## 
## [[48]]
## [1] 1.547248e-04 6.415608e-05 9.502806e-05 5.052029e-04
## 
## [[49]]
##  [1] 1.165162e-04 1.276375e-04 8.762809e-05 6.065400e-05 6.157729e-05
##  [6] 2.138680e-04 6.349008e-05 8.240884e-05 1.846241e-04 6.681049e-05
## [11] 7.612655e-05
## 
## [[50]]
## [1] 1.016815e-04 7.644401e-05 9.641867e-05 1.096279e-04
## 
## [[51]]
## [1] 1.161859e-04 2.417046e-04 9.240346e-05 3.270855e-04
## 
## [[52]]
## [1] 1.728270e-04 7.494096e-05 1.923145e-04
## 
## [[53]]
## [1] 1.741204e-04 5.500956e-05 6.797303e-05 9.671140e-05 7.848568e-05
## 
## [[54]]
## [1] 1.092489e-04 8.780628e-05 9.145785e-05 1.069529e-04 1.119878e-04
## 
## [[55]]
## [1] 1.414978e-04 1.292529e-04 5.954553e-05 1.312002e-04 2.425602e-04
## [6] 6.763797e-05
## 
## [[56]]
## [1] 6.851713e-05 1.426785e-04 1.322068e-04 4.494805e-04 1.397649e-04
## 
## [[57]]
## [1] 6.415608e-05 1.996281e-04 8.649248e-05
## 
## [[58]]
## [1] 1.274074e-04 7.998431e-05 1.014577e-04 7.556778e-05 2.063834e-04
## 
## [[59]]
## [1] 2.639822e-04 4.687696e-04 7.735722e-05 6.138432e-05
## 
## [[60]]
## [1] 6.157729e-05 5.644105e-05 1.253675e-04 1.859244e-04 1.274989e-04
## [6] 1.228858e-04
## 
## [[61]]
## [1] 5.074741e-05 2.958042e-04 9.027376e-05 6.813486e-05 8.333363e-05
## 
## [[62]]
## [1] 1.355575e-04 8.780628e-05 7.087537e-05 1.078882e-04 6.211540e-05
## 
## [[63]]
## [1] 7.998431e-05 5.845233e-05 1.097591e-04 1.028617e-04 5.773333e-05
## [6] 1.110776e-04
## 
## [[64]]
## [1] 9.289630e-05 1.014577e-04 5.845233e-05 7.342780e-05 6.392337e-05
## [6] 8.015172e-05 6.542939e-05
## 
## [[65]]
## [1] 6.064992e-05 2.347901e-04 3.899482e-04 1.669993e-04 9.817634e-05
## [6] 7.966890e-05
## 
## [[66]]
## [1] 1.178655e-04 1.126621e-04 5.261504e-05 5.840214e-05 6.909403e-05
## [6] 8.609568e-05 1.094156e-04 1.008985e-04 6.612251e-05
## 
## [[67]]
##  [1] 6.333917e-05 6.797303e-05 5.479950e-05 6.819981e-05 2.501014e-04
##  [6] 3.845427e-04 8.021070e-05 1.622446e-04 4.930344e-04 1.382934e-04
## 
## [[68]]
## [1] 6.064992e-05 7.694168e-05 6.040944e-05 7.272001e-05 4.471832e-04
## [6] 1.427255e-04
## 
## [[69]]
## [1] 6.310945e-05 5.713328e-05 1.525564e-04
## 
## [[70]]
## [1] 5.982625e-05 6.367085e-05 5.479950e-05 1.653097e-04 1.252556e-04
## [6] 8.165058e-05 5.846018e-05
## 
## [[71]]
## [1] 6.981223e-05 6.277917e-05 1.653097e-04
## 
## [[72]]
## [1] 1.570112e-04 7.717293e-05 1.097591e-04 7.342780e-05 1.147446e-04
## [6] 7.717540e-05
## 
## [[73]]
## [1] 3.382463e-04 8.881606e-05 1.246360e-04 9.671140e-05 7.087537e-05
## [6] 6.819981e-05 6.120162e-05 1.023066e-04
## 
## [[74]]
## [1] 2.448648e-04 9.734377e-05 1.447864e-04 5.586925e-05 5.328770e-05
## [6] 2.138680e-04 1.676493e-04 1.235205e-04
## 
## [[75]]
## [1] 1.502051e-04 5.644105e-05 5.713328e-05 7.671582e-05 8.538918e-05
## [6] 5.975063e-05 7.195585e-05 7.104666e-05
## 
## [[76]]
## [1] 6.797692e-05 9.145785e-05 1.028617e-04 8.037875e-05 8.662185e-05
## [6] 8.293435e-05 1.404467e-04
## 
## [[77]]
## [1] 0.0001918575
## 
## [[78]]
## [1] 6.317576e-05 7.378774e-05 6.829824e-05 5.920533e-05 8.473366e-05
## [6] 7.848568e-05 2.501014e-04 4.617344e-04
## 
## [[79]]
## [1] 4.687696e-04 2.347901e-04 7.694168e-05 1.490246e-04 7.638188e-05
## [6] 3.000689e-04 6.315298e-05 8.856767e-05
## 
## [[80]]
## [1] 6.392337e-05 8.794944e-05 1.108574e-04 6.711318e-05
## 
## [[81]]
## [1] 6.965538e-05 9.734742e-05 2.243529e-04 1.687366e-04
## 
## [[82]]
## [1] 1.529304e-04 1.173296e-04 8.839082e-05 4.549339e-04 6.245430e-05
## 
## [[83]]
## [1] 1.069529e-04 1.078882e-04 3.845427e-04 5.643471e-05 2.420899e-04
## [6] 1.193623e-04
## 
## [[84]]
## [1] 7.456176e-05 5.614833e-05 1.142536e-04 5.840214e-05 1.178126e-04
## [6] 1.427576e-04 9.904001e-05
## 
## [[85]]
## [1] 7.616479e-05 7.671596e-05 2.505499e-04
## 
## [[86]]
## [1] 7.795105e-05 1.071633e-04 7.801567e-05 5.954553e-05 1.676493e-04
## [6] 4.374973e-04
## 
## [[87]]
## [1] 1.147090e-04 9.027376e-05
## 
## [[88]]
## [1] 7.253707e-05 1.956579e-04 1.990647e-04 1.073083e-04
## 
## [[89]]
## [1] 7.591507e-05 8.074777e-05 6.435784e-05 1.956579e-04 8.824106e-05
## [6] 7.768063e-05
## 
## [[90]]
## [1] 1.794158e-04 7.556778e-05 8.015172e-05 8.794944e-05 5.870988e-05
## [6] 7.304769e-05
## 
## [[91]]
## [1] 1.148537e-04 9.734742e-05 7.392934e-05 9.734426e-05 6.873991e-05
## 
## [[92]]
## [1] 8.021070e-05 1.252556e-04 6.898998e-05
## 
## [[93]]
## [1] 9.521408e-05 1.103750e-04 2.486247e-04
## 
## [[94]]
## [1] 0.0001996281 0.0001529304 0.0001467189
## 
## [[95]]
## [1] 6.813486e-05
## 
## [[96]]
## [1] 1.989998e-04 9.839306e-05 1.654548e-04 1.710716e-04 1.119231e-04
## 
## [[97]]
## [1] 9.455630e-05 6.349008e-05 1.253675e-04 7.671582e-05 8.285562e-05
## [6] 8.432010e-05 9.803931e-05 1.269201e-04
## 
## [[98]]
## [1] 6.215702e-05 9.521408e-05 1.654548e-04 1.534867e-04 9.459651e-05
## [6] 9.267810e-05 9.259669e-05
## 
## [[99]]
## [1] 5.348017e-05 1.544491e-04 6.909403e-05 7.184076e-05 2.280988e-04
## 
## [[100]]
## [1] 1.461145e-04 7.392934e-05 8.285562e-05 6.810813e-05
## 
## [[101]]
## [1] 5.669883e-05 7.459317e-05 6.542939e-05 1.108574e-04
## 
## [[102]]
## [1] 7.573288e-05 1.162834e-04 1.173296e-04 1.103750e-04 1.534867e-04
## [6] 1.684752e-04 2.260939e-04 1.430287e-04 9.710852e-05
## 
## [[103]]
## [1] 9.416121e-05 1.119878e-04 2.063834e-04 5.773333e-05 8.037875e-05
## [6] 5.643471e-05 5.936774e-05
## 
## [[104]]
## [1] 3.899482e-04 6.040944e-05 1.490246e-04 2.629296e-04 6.504732e-05
## [6] 2.120856e-04
## 
## [[105]]
## [1] 2.320960e-04 1.341811e-04 8.609568e-05 7.184076e-05 6.360556e-05
## 
## [[106]]
## [1] 8.143763e-05 8.780240e-05 1.669993e-04 1.990647e-04 2.629296e-04
## [6] 7.205777e-05 7.535666e-05 5.727673e-05 2.571497e-04
## 
## [[107]]
## [1] 0.0001205144 0.0000504364 0.0002379534 0.0003872600
## 
## [[108]]
## [1] 6.162470e-05 2.243529e-04 9.734426e-05 8.432010e-05 6.810813e-05
## [6] 2.338874e-04
## 
## [[109]]
## [1] 7.825048e-05 1.400602e-04 1.687366e-04 7.671596e-05 6.873991e-05
## [6] 2.338874e-04 6.080557e-05
## 
## [[110]]
## [1] 3.757265e-03 1.094156e-04 1.178126e-04 8.142763e-05 6.656400e-05
## 
## [[111]]
## [1] 8.003154e-05 6.274283e-05 2.417046e-04 8.662185e-05 2.379534e-04
## [6] 1.299429e-04 1.105673e-04
## 
## [[112]]
## [1] 5.240886e-05 9.240346e-05 7.494096e-05 1.227695e-04
## 
## [[113]]
## [1] 2.572725e-04 6.920613e-05 2.486247e-04 1.710716e-04 9.459651e-05
## [6] 1.684752e-04
## 
## [[114]]
## [1] 2.928491e-04 9.817634e-05 1.073083e-04 8.824106e-05 7.205777e-05
## [6] 3.641001e-04
## 
## [[115]]
##  [1] 1.489506e-04 8.314604e-05 6.943647e-05 5.725829e-05 7.605560e-05
##  [6] 8.240884e-05 1.312002e-04 1.235205e-04 2.505499e-04 6.080557e-05
## 
## [[116]]
## [1] 7.484700e-05 6.348349e-05 6.054047e-04
## 
## [[117]]
##  [1] 6.716088e-05 2.055733e-04 1.046676e-04 9.297839e-05 1.036850e-04
##  [6] 9.400085e-05 7.735722e-05 7.966890e-05 1.622446e-04 8.165058e-05
## [11] 4.617344e-04 7.638188e-05 7.768063e-05 3.641001e-04
## 
## [[118]]
## [1] 1.063840e-04 1.232539e-04 6.827966e-05 1.130843e-04 8.333363e-05
## [6] 1.147446e-04
## 
## [[119]]
## [1] 7.246207e-05 5.013834e-05 5.986237e-05 1.053480e-04 2.425602e-04
## [6] 2.280988e-04 6.360556e-05
## 
## [[120]]
## [1] 7.435495e-05 7.937698e-05 1.155047e-04 1.846241e-04 6.120162e-05
## [6] 9.217043e-05 6.066269e-05
## 
## [[121]]
## [1] 7.026734e-05
## 
## [[122]]
## [1] 1.114368e-04 1.617230e-04 4.930344e-04 5.846018e-05 2.420899e-04
## [6] 6.898998e-05 5.936774e-05 7.535666e-05 6.060526e-05
## 
## [[123]]
## [1] 0.0001426785 0.0000709919
## 
## [[124]]
## [1] 7.768842e-05 6.583885e-05 1.525564e-04 8.538918e-05 1.768949e-04
## 
## [[125]]
##  [1] 9.351276e-05 1.265698e-04 6.063476e-05 6.681049e-05 1.859244e-04
##  [6] 5.975063e-05 9.217043e-05 1.768949e-04 8.113619e-05 3.582431e-04
## 
## [[126]]
## [1] 6.228594e-05 9.502806e-05 8.649248e-05 8.839082e-05 1.467189e-04
## [6] 2.260939e-04 9.371913e-05
## 
## [[127]]
## [1] 1.206285e-04 5.052029e-04 5.969137e-05 6.386635e-05
## 
## [[128]]
## [1] 1.658591e-04 1.110776e-04 7.717540e-05 8.293435e-05 1.299429e-04
## 
## [[129]]
## [1] 1.108079e-04 1.322068e-04 6.504732e-05 5.727673e-05 6.794786e-05
## 
## [[130]]
## [1] 7.500310e-05 1.150039e-04 5.870988e-05 5.969137e-05 2.783276e-04
## 
## [[131]]
## [1] 0.0000621154 0.0001382934 0.0001023066 0.0001193623
## 
## [[132]]
## [1] 4.549339e-04 9.267810e-05 1.430287e-04 9.371913e-05 5.634864e-05
## 
## [[133]]
## [1] 1.955575e-04 6.175376e-05 4.494805e-04 1.119231e-04
## 
## [[134]]
## [1] 1.043238e-04 8.822855e-05 6.794225e-05 9.710852e-05 2.571497e-04
## [6] 6.060526e-05
## 
## [[135]]
## [1] 1.397649e-04 7.272001e-05 2.120856e-04 7.099190e-05 6.794786e-05
## [6] 6.353800e-05
## 
## [[136]]
## [1] 1.749139e-04 5.785814e-05 3.270855e-04 1.923145e-04 1.274989e-04
## [6] 7.195585e-05 9.803931e-05 3.872600e-04 1.227695e-04
## 
## [[137]]
## [1] 6.979176e-05 1.370508e-04 1.404467e-04 1.105673e-04 6.066269e-05
## [6] 8.113619e-05
## 
## [[138]]
## [1] 7.082767e-05 1.367894e-04 9.422543e-05 6.763797e-05 4.374973e-04
## 
## [[139]]
## [1] 6.245430e-05 9.259669e-05 5.634864e-05
## 
## [[140]]
## [1] 6.711318e-05 7.304769e-05 6.386635e-05 2.783276e-04
## 
## [[141]]
## [1] 7.612655e-05 1.228858e-04 7.104666e-05 1.269201e-04 3.582431e-04
## 
## [[142]]
## [1] 1.241413e-04 8.539523e-05 6.672324e-05 6.983264e-05
## 
## [[143]]
## [1] 1.427576e-04 8.142763e-05 5.544160e-05 5.829138e-05 3.362260e-04
## [6] 1.781684e-04 1.924510e-04
## 
## [[144]]
## [1] 1.243936e-04 9.994883e-05 9.129165e-05
## 
## [[145]]
## [1] 3.281936e-04 1.928098e-04 5.898352e-04 3.083782e-04 7.644762e-05
## [6] 7.580708e-05 9.554232e-05 6.357697e-05 7.440387e-05
## 
## [[146]]
## [1] 4.471832e-04 6.353800e-05 1.241413e-04 6.544388e-05 1.824225e-04
## [6] 1.365450e-04
## 
## [[147]]
## [1] 5.544160e-05 3.281936e-04 7.409936e-05 6.044252e-05
## 
## [[148]]
## [1] 3.000689e-04 1.656710e-04 5.674852e-05 1.314863e-04 8.146334e-05
## [6] 6.177564e-05 3.450862e-04
## 
## [[149]]
## [1] 1.656710e-04 1.640892e-04 7.444357e-05 2.661668e-04 3.303072e-04
## [6] 1.041609e-04
## 
## [[150]]
## [1] 6.381720e-05 7.704263e-05 6.138432e-05 6.315298e-05 5.674852e-05
## [6] 6.730039e-05
## 
## [[151]]
## [1] 1.243936e-04 1.928098e-04 2.863397e-04 5.493352e-05
## 
## [[152]]
## [1] 0.0001314863 0.0001640892 0.0001973662 0.0002709325
## 
## [[153]]
## [1] 7.283566e-05 8.508572e-05 2.643966e-04 8.697488e-05
## 
## [[154]]
## [1] 9.994883e-05 5.898352e-04 2.863397e-04 2.126969e-04
## 
## [[155]]
## [1] 8.539523e-05 6.544388e-05 9.720004e-05 2.421669e-04 1.030413e-04
## [6] 1.613814e-04
## 
## [[156]]
## [1] 1.427255e-04 8.856767e-05 1.824225e-04 8.146334e-05 9.720004e-05
## [6] 1.918939e-04 1.342323e-04
## 
## [[157]]
## [1] 7.644401e-05 1.385325e-04 7.110111e-05 2.916675e-04
## 
## [[158]]
## [1] 6.672324e-05 7.283566e-05 2.421669e-04 7.108129e-05 1.008078e-04
## [6] 6.897220e-05
## 
## [[159]]
## [1] 6.177564e-05 7.444357e-05 1.918939e-04 6.121227e-05 1.094021e-04
## 
## [[160]]
## [1] 9.129165e-05 3.083782e-04 2.126969e-04 9.969477e-05
## 
## [[161]]
## [1] 1.259769e-04 3.450862e-04 6.730039e-05 1.973662e-04 8.760308e-05
## [6] 6.111340e-05 1.766843e-04
## 
## [[162]]
## [1] 6.983264e-05 1.365450e-04 1.030413e-04
## 
## [[163]]
## [1] 3.042400e-04 9.641867e-05 1.385325e-04 1.668645e-04
## 
## [[164]]
## [1] 1.883358e-04 1.008985e-04 8.760308e-05 2.016158e-04 9.397873e-05
## 
## [[165]]
## [1] 7.644762e-05 2.661668e-04 2.636170e-04 2.315066e-04 1.279681e-04
## [6] 8.041173e-05 2.298781e-04
## 
## [[166]]
## [1] 7.580708e-05 8.508572e-05 7.108129e-05 9.969477e-05 2.636170e-04
## [6] 9.518302e-05
## 
## [[167]]
## [1] 5.829138e-05 6.111340e-05 2.016158e-04 6.464456e-05 5.803879e-05
## 
## [[168]]
## [1] 6.612251e-05 6.656400e-05 3.362260e-04 9.397873e-05 6.464456e-05
## 
## [[169]]
## [1] 3.303072e-04 6.121227e-05 2.315066e-04 1.018157e-04 1.452857e-04
## 
## [[170]]
## [1] 1.781684e-04 9.554232e-05 7.409936e-05 1.041609e-04 2.709325e-04
## [6] 1.766843e-04 1.279681e-04 5.803879e-05
## 
## [[171]]
## [1] 9.680763e-05 1.096279e-04 9.904001e-05 1.924510e-04 6.357697e-05
## [6] 6.044252e-05 7.110111e-05 8.893700e-05 1.067492e-04
## 
## [[172]]
## [1] 2.643966e-04 8.041173e-05 9.518302e-05 9.999673e-05
## 
## [[173]]
## [1] 0.0002916675 0.0001668645 0.0000889370 0.0000619717
## 
## [[174]]
## [1] 7.440387e-05 5.493352e-05 1.067492e-04 6.197170e-05
## 
## [[175]]
## [1] 0.0001613814 0.0001342323 0.0001008078 0.0001094021 0.0001018157
## [6] 0.0001244095
## 
## [[176]]
## [1] 8.697488e-05 6.897220e-05 2.298781e-04 1.452857e-04 9.999673e-05
## [6] 1.244095e-04

6.3 Row-standardised weights matrix

Next, we need to assign weights to each neighboring polygon. In our case, each neighboring polygon will be assigned equal weight (style=“W”). This is accomplished by assigning the fraction 1/(#ofneighbors) to each neighboring county then summing the weighted income values. While this is the most intuitive way to summaries the neighbors’ values it has one drawback in that polygons along the edges of the study area will base their lagged values on fewer polygons thus potentially over- or under-estimating the true nature of the spatial autocorrelation in the data. For this example, we’ll stick with the style=“W” option for simplicity’s sake but note that other more robust options are available, notably style=“B”.

rswm_q_mexi_combined_filtered <- nb2listw(wm_q_mexi_combined_filtered, style="W", zero.policy = TRUE)
rswm_q_mexi_combined_filtered
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 962 
## Percentage nonzero weights: 3.10563 
## Average number of links: 5.465909 
## 
## Weights style: W 
## Weights constants summary:
##     n    nn  S0       S1      S2
## W 176 30976 176 71.10695 731.679
rswm_ids <- nb2listw(wm_q_mexi_combined_filtered, glist=ids, style="B", zero.policy=TRUE)
rswm_ids
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 962 
## Percentage nonzero weights: 3.10563 
## Average number of links: 5.465909 
## 
## Weights style: B 
## Weights constants summary:
##     n    nn        S0           S1           S2
## B 176 30976 0.1314899 0.0001067155 0.0005739984
summary(unlist(rswm_ids$weights))
##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## 5.014e-05 7.083e-05 9.422e-05 1.367e-04 1.461e-04 3.757e-03

7 Application of Spatial Weight Matrix

Using spatial lag variables with row-standardized weights,

COVID19RATE.lag <- lag.listw(rswm_q_mexi_combined_filtered, mexi_combined_filtered$rate_of_covid)
COVID19RATE.lag
##   [1]  65.173501  83.919003  56.334680  58.841201  69.311153  69.831564
##   [7]  66.897849  54.460995  78.951095  86.628452  72.977175 105.188880
##  [13]  80.881738  86.497243  71.390864  74.573109  18.430921  27.737866
##  [19]  14.423635   8.818124  24.422779  60.137774  16.501641  13.762601
##  [25]  23.508170  19.893033  27.552122  35.310943  27.892150  16.669582
##  [31]  20.322240  38.075314  25.812475  41.208353  29.509486  44.907087
##  [37]  17.949159  54.472268  25.413166  27.721732  57.078502  14.590662
##  [43]  40.425085  22.993079  26.799172  21.897430  29.690661  19.024676
##  [49]  34.577960  14.085890  25.902176  23.436928  56.906521  22.877109
##  [55]  31.345833  12.464604  19.622254  16.927449  37.071844  30.400805
##  [61]  11.155867  29.270920  18.313267  20.654229  24.118175  45.444376
##  [67]  36.069927  21.110282  38.100668  40.429342  42.517254  18.208294
##  [73]  49.650651  56.415698  34.063709  16.019590  40.945406  62.080004
##  [79]  21.346132  15.835091  32.250376  19.406752  30.280135  17.260040
##  [85]  18.893468  34.496693  13.015758  31.259766  54.609969  18.707917
##  [91]  31.262303  41.814313  20.893620  15.992568  19.706406  16.481862
##  [97]  30.716693  14.663874  35.155219  37.274139  28.293344  22.702941
## [103]  21.290845  18.542540  37.226166  27.497065  36.935359  26.366913
## [109]  29.803270  14.772992  29.605548  24.473185  17.617340  32.969111
## [115]  25.948673  22.300648  46.131139  16.854869  38.624164  49.540459
## [121]   3.519392  29.252209  18.770633  36.949700  37.596190  13.772708
## [127]  14.287846  14.442223  21.209975  15.718976  30.525001  22.537199
## [133]  12.476061  23.599802  18.210045  29.324160  32.713512  58.209694
## [139]  19.037816  17.479400  28.616720  20.669144  22.456103  14.732444
## [145]  16.416916  14.048310  15.693980  14.703250  17.335272  57.362738
## [151]  17.793050  15.452228  25.844834  19.594349  15.897136  14.668824
## [157]  12.370439  12.968487  14.771576  21.108920  42.558303  13.548369
## [163]  12.081496  48.242039  18.675616  27.655608  13.932382  19.863409
## [169]  10.830841  20.794401  16.323098  18.015965  11.669421  13.576173
## [175]  17.352059  23.677386

Append the spatially lag COVID19RATE.lag values onto mexi_combined_filtered SpatialPolygonDataFrame by using the code chunk below.

lag.list <- list(mexi_combined_filtered$NOMGEO, lag.listw(rswm_q_mexi_combined_filtered, mexi_combined_filtered$rate_of_covid))
lag.res <- as.data.frame(lag.list)
colnames(lag.res) <- c("NOMGEO", "lag COVID 19 Rate")
mexi_combined_filtered@data <- left_join(mexi_combined_filtered@data,lag.res)
head(mexi_combined_filtered@data, 10)
##    CVEGEO CVE_ENT CVE_MUN                 NOMGEO Pop2010 Pop2020 new13 new14
## 1   09002      09     002           Azcapotzalco  414711  408441    18    27
## 2   09003      09     003              Coyoacán  620416  621952    15    19
## 3   09004      09     004  Cuajimalpa de Morelos  186391  199809    18     9
## 4   09005      09     005      Gustavo A. Madero 1185772 1176967    21    68
## 5   09006      09     006              Iztacalco  384326  393821     9    16
## 6   09007      09     007             Iztapalapa 1815786 1815551    32    59
## 7   09008      09     008 La Magdalena Contreras  239086  245147     9    12
## 8   09009      09     009             Milpa Alta  130582  139371     0     1
## 9   09010      09     010       Ã\201lvaro Obregón  727034  755537    18    34
## 10  09011      09     011               Tláhuac  360265  366586     2     8
##    new15 new16 new17 new18 new19 new20 new21 new22 new23 new24 new25 new26
## 1     41    69    98   137   224   316   336   323   349   310   345   328
## 2     63   106   162   179   277   319   342   374   389   299   382   272
## 3     18    29    35    45    58    63    91   159   141   137   140   126
## 4     97   204   337   373   546   677   778   735   831   692   611   525
## 5     38   123   155   209   254   319   305   276   283   220   245   206
## 6    125   301   535   693   863   966  1076   873   884   818   747   665
## 7     18    38    86    51    74   133   144   142   126   122   181   179
## 8     22    36    78    70    94   119   195   163   183   223   170   137
## 9     56   124   142   197   238   322   385   387   383   458   435   367
## 10    19    49    96   135   203   303   352   233   303   237   279   277
##    new27 new28 new29 new30 new31 new32 cumul13 cumul14 cumul15 cumul16 cumul17
## 1    316   301   300   332   271    24      19      46      87     156     254
## 2    382   313   490   427   396    30      25      44     107     213     375
## 3    168   159   171   220   114     7      66      75      93     122     157
## 4    573   512   625   663   505    59      28      96     193     397     734
## 5    165   172   186   233   196    62      11      27      65     188     343
## 6    762   670   744   728   620   135      41     100     225     526    1061
## 7    208   229   345   247   277    12      15      27      45      83     169
## 8    176   151   204   197   145    11       0       1      23      59     137
## 9    388   419   564   492   425    62      56      90     146     270     412
## 10   295   312   351   326   285    38       3      11      30      79     175
##    cumul18 cumul19 cumul20 cumul21 cumul22 cumul23 cumul24 cumul25 cumul26
## 1      391     615     931    1267    1590    1939    2249    2594    2922
## 2      554     831    1150    1492    1866    2255    2554    2936    3208
## 3      202     260     323     414     573     714     851     991    1117
## 4     1107    1653    2330    3108    3843    4674    5366    5977    6502
## 5      552     806    1125    1430    1706    1989    2209    2454    2660
## 6     1754    2617    3583    4659    5532    6416    7234    7981    8646
## 7      220     294     427     571     713     839     961    1142    1321
## 8      207     301     420     615     778     961    1184    1354    1491
## 9      609     847    1169    1554    1941    2324    2782    3217    3584
## 10     310     513     816    1168    1401    1704    1941    2220    2497
##    cumul27 cumul28 cumul29 cumul30 cumul31 cumul32 activ13 activ14 activ15
## 1     3238    3539    3839    4171    4442    4466      34      70     130
## 2     3590    3903    4393    4820    5216    5246      41      82     166
## 3     1285    1444    1615    1835    1949    1956      70      86     112
## 4     7075    7587    8212    8875    9380    9439      68     173     321
## 5     2825    2997    3183    3416    3612    3674      20      58     129
## 6     9408   10078   10822   11550   12170   12305      75     174     427
## 7     1529    1758    2103    2350    2627    2639      24      41      70
## 8     1667    1818    2022    2219    2364    2375       1       9      43
## 9     3972    4391    4955    5447    5872    5934      80     135     215
## 10    2792    3104    3455    3781    4066    4104       8      29      61
##    activ16 activ17 activ18 activ19 activ20 activ21 activ22 activ23 activ24
## 1      199     332     516     776    1126    1456    1795    2143    2451
## 2      321     513     731    1054    1398    1751    2117    2496    2841
## 3      140     182     237     298     385     498     660     801     917
## 4      571     940    1438    2015    2755    3557    4275    5115    5756
## 5      286     482     722    1008    1328    1618    1893    2147    2394
## 6      854    1492    2315    3207    4207    5208    6095    6979    7723
## 7      113     206     280     386     507     653     786     925    1057
## 8      107     186     288     388     548     737     921    1129    1304
## 9      368     545     758    1053    1403    1807    2172    2651    3058
## 10     142     257     441     695    1023    1317    1579    1881    2111
##    activ25 activ26 activ27 activ28 activ29 activ30 activ31 activ32 death13
## 1     2810    3117    3444    3746    4034    4340    4460    4466       0
## 2     3169    3496    3862    4227    4681    5075    5245    5246       0
## 3     1080    1181    1384    1540    1677    1885    1954    1956       0
## 4     6329    6842    7459    8022    8636    9182    9427    9439       2
## 5     2591    2766    2953    3137    3369    3553    3666    3674       1
## 6     8434    9100    9874   10528   11287   11950   12283   12305       2
## 7     1221    1408    1626    1877    2170    2446    2638    2639       0
## 8     1481    1626    1787    1964    2170    2322    2375    2375       0
## 9     3471    3812    4286    4744    5284    5705    5924    5934       2
## 10    2408    2681    3016    3325    3664    3940    4101    4104       0
##    death14 death15 death16 death17 death18 death19 death20 death21 death22
## 1        2       4       6      15      14      38      33      44      47
## 2        1       7       8      19      11      26      30      36      34
## 3        2       1       3       5       1       4      11       8       9
## 4        5      26      19      37      74      82      98     117     114
## 5        1       3       7      26      20      30      46      49      37
## 6        6      15      33      76     115     127     174     176     127
## 7        1       4       3       0       9       4      16      12       5
## 8        0       0       1       3       5       7       4       3       7
## 9        2      11      11      23      29      31      50      40      60
## 10       1       1       2       5      16      14      14      22      17
##    death23 death24 death25 death26 death27 death28 death29 death30 death31
## 1       44      43      47      38      27      19      31      20      19
## 2       41      34      37      33      31      25      23      20      16
## 3       11      18       6       5       5       9      10       7       8
## 4      121     148     112      77      72      56      57      52      33
## 5       29      35      26      22      19      25      23      17      12
## 6      120     110     103      63      67      65      52      29      29
## 7       19      10      10       7      12       8       5       5       6
## 8        4       6       7      11       4       2       2       0       2
## 9       56      61      41      32      39      26      28      20      24
## 10      20      17      24      12      11       8      10       3       8
##    death32 wk13_32_cumul_count wk13_32_activ_count wk13_32_death_count
## 1        8                4465               41445                 499
## 2        6                5236               48512                 438
## 3        1                1908               17043                 124
## 4       18                9432               92320                1320
## 5        3                3672               37794                 431
## 6       17               12296              124517                1506
## 7        3                2633               21073                 139
## 8        0                2375               21761                  68
## 9        8                5896               53405                 594
## 10       2                4103               36783                 207
##    rate_of_covid active_rate_of_covid death_rate_of_covid lag COVID 19 Rate
## 1      109.31811            1014.7120           12.217187          65.17350
## 2       84.18656             779.9959            7.042344          83.91900
## 3       95.49119             852.9646            6.205927          56.33468
## 4       80.13819             784.3890           11.215268          58.84120
## 5       93.24033             959.6746           10.944058          69.31115
## 6       67.72600             685.8359            8.295002          69.83156
## 7      107.40494             859.6067            5.670067          66.89785
## 8      170.40848            1561.3722            4.879064          54.46100
## 9       78.03721             706.8482            7.861958          78.95109
## 10     111.92462            1003.3935            5.646697          86.62845

Plotting both the COVID19Rate and spatial lag COVID19Rate for comparison using the code chunk below.

covid <- qtm(mexi_combined_filtered, "rate_of_covid")
lag_covid <- qtm(mexi_combined_filtered, "lag COVID 19 Rate")
tmap_arrange(covid, lag_covid, asp=1, ncol=2)

# Measure of Global Spatial Autocorrelation

7.1 Maron’s I test

Perform Moran’s I statistics testing using moran.test() of spdep.

moran.test(mexi_combined_filtered$rate_of_covid, listw=rswm_q_mexi_combined_filtered, zero.policy = TRUE, na.action=na.omit)
## 
##  Moran I test under randomisation
## 
## data:  mexi_combined_filtered$rate_of_covid  
## weights: rswm_q_mexi_combined_filtered    
## 
## Moran I statistic standard deviate = 10.854, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.495940195      -0.005714286       0.002135988

The Moran I statistic shows that the value of 0.495940195 is a positive spatial autocorrelation that of a sign of Clustering. The p value we got is lesser than the targeted 0.05 does we reject the null hypothesis where distribution of COVID 19 rate is randomly distributed.

7.2 Computing Monte Carlo Moran’s I

The code chunk below performs permutation test for Moran’s I statistic by using moran.mc() of spdep.

set.seed(1234)
bperm= moran.mc(mexi_combined_filtered$rate_of_covid, listw=rswm_q_mexi_combined_filtered, nsim=999, zero.policy = TRUE, na.action=na.omit)
bperm
## 
##  Monte-Carlo simulation of Moran I
## 
## data:  mexi_combined_filtered$rate_of_covid 
## weights: rswm_q_mexi_combined_filtered  
## number of simulations + 1: 1000 
## 
## statistic = 0.49594, observed rank = 1000, p-value = 0.001
## alternative hypothesis: greater

7.3 Visualising Monte Carlo Moran’s I

mean(bperm$res[1:999])
## [1] -0.006244441
var(bperm$res[1:999])
## [1] 0.002138681
summary(bperm$res[1:999])
##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## -0.133572 -0.039006 -0.009037 -0.006244  0.023391  0.193374
hist(bperm$res, freq=TRUE, breaks=50, xlab="Simulated Moran's I")
abline(v=0, col="red") 

## Computing local Moran’s I

To compute local Moran’s I, the localmoran() function of spdep will be used. It computes Ii values, given a set of zi values and a listw object providing neighbour weighting information for the polygon associated with the zi values.

fips <- order(mexi_combined_filtered$NOMGEO)
localMI <- localmoran(mexi_combined_filtered$rate_of_covid, rswm_q_mexi_combined_filtered)
head(localMI)
##           Ii         E.Ii    Var.Ii     Z.Ii    Pr(z > 0)
## 270 3.715790 -0.005714286 0.1857254 8.635408 2.925889e-18
## 271 3.907305 -0.005714286 0.1857254 9.079802 5.438625e-20
## 272 2.268934 -0.005714286 0.2334478 4.707815 1.251934e-06
## 273 1.899870 -0.005714286 0.1141417 5.640347 8.485396e-09
## 274 3.310531 -0.005714286 0.1857254 7.695043 7.072335e-15
## 275 1.980378 -0.005714286 0.1141417 5.878644 2.068199e-09

localmoran() function returns a matrix of values whose columns are:

Ii: the local Moran’s I statistics E.Ii: the expectation of local moran statistic under the randomisation hypothesis Var.Ii: the variance of local moran statistic under the randomisation hypothesis Z.Ii:the standard deviate of local moran statistic Pr(): the p-value of local moran statistic The code chunk below list the content of the local Moran matrix derived by using printCoefmat().

printCoefmat(data.frame(localMI[fips,], row.names=mexi_combined_filtered$NOMGEO[fips]), check.names=FALSE)
##                                        Ii        E.Ii      Var.Ii        Z.Ii
## Ã\201lvaro Obregón               3.13141221 -0.00571429  0.13118547  8.66142143
## Acambay de Ruíz Castañeda    0.26390635 -0.00571429  0.23344785  0.55803082
## Acolman                        0.01426135 -0.00571429  0.13118547  0.05515157
## Aculco                         0.56632612 -0.00571429  0.23344785  1.18394565
## Almoloya de Alquisiras         0.30401140 -0.00571429  0.23344785  0.64103579
## Almoloya de Juárez            0.13581163 -0.00571429  0.13118547  0.39074473
## Almoloya del Río              0.67470812 -0.00571429  0.23344785  1.40826267
## Amacuzac                       0.32183542 -0.00571429  0.23344785  0.67792598
## Amanalco                      -0.01484814 -0.00571429  0.13118547 -0.02521803
## Amatepec                       0.65595393 -0.00571429  0.31298526  1.18271017
## Amecameca                     -0.05132161 -0.00571429  0.23344785 -0.09439296
## Apaxco                         0.19887632 -0.00571429  0.47206007  0.29777420
## Atenco                         0.05565255 -0.00571429  0.15391044  0.15642267
## Atizapán                      0.74051905 -0.00571429  0.31298526  1.33386755
## Atizapán de Zaragoza          0.02399602 -0.00571429  0.15391044  0.07573090
## Atlacomulco                   -0.11524758 -0.00571429  0.18572540 -0.25416192
## Atlatlahucan                   0.21459988 -0.00571429  0.13118547  0.60827443
## Atlautla                       0.24540192 -0.00571429  0.18572540  0.58269205
## Axapusco                       0.09806448 -0.00571429  0.23344785  0.21478976
## Axochiapan                     0.05576287 -0.00571429  0.31298526  0.10988839
## Ayala                          0.18529881 -0.00571429  0.10088550  0.60137970
## Ayapango                       0.12507698 -0.00571429  0.15391044  0.33338396
## Azcapotzalco                   3.71578995 -0.00571429  0.18572540  8.63540803
## Benito Juárez                 2.04673506 -0.00571429  0.15391044  5.23164671
## Calimaya                       0.02624586 -0.00571429  0.11414173  0.09459899
## Capulhuac                     -0.00302759 -0.00571429  0.31298526  0.00480237
## Chalco                         0.92752879 -0.00571429  0.10088550  2.93819355
## Chapa de Mota                  0.46452020 -0.00571429  0.23344785  0.97323907
## Chapultepec                    0.23978916 -0.00571429  0.23344785  0.50811575
## Chiautla                       0.12379509 -0.00571429  0.13118547  0.35756776
## Chicoloapan                    0.05317621 -0.00571429  0.23344785  0.12188501
## Chiconcuac                     0.09787837 -0.00571429  0.31298526  0.18516847
## Chimalhuacán                  0.00927770 -0.00571429  0.23344785  0.03102874
## Coacalco de Berriozábal       0.04804869 -0.00571429  0.23344785  0.11127263
## Coatepec Harinas               0.45174810 -0.00571429  0.11414173  1.35404518
## Coatetelco                     0.55220367 -0.00571429  0.15391044  1.42212017
## Coatlán del Río              0.30183344 -0.00571429  0.15391044  0.78393217
## Cocotitlán                    0.43055795 -0.00571429  0.31298526  0.77982227
## Coyoacán                      3.90730541 -0.00571429  0.18572540  9.07980204
## Coyotepec                     -0.04554039 -0.00571429  0.23344785 -0.08242763
## Cuajimalpa de Morelos          2.26893414 -0.00571429  0.23344785  4.70781452
## Cuauhtémoc                    3.75433271 -0.00571429  0.15391044  9.58427430
## Cuautitlán                   -0.13585213 -0.00571429  0.11414173 -0.38519565
## Cuautitlán Izcalli           -0.02486671 -0.00571429  0.15391044 -0.04881909
## Cuautla                       -0.06787939 -0.00571429  0.23344785 -0.12866242
## Cuernavaca                     0.11106901 -0.00571429  0.13118547  0.32243180
## Donato Guerra                  0.37624778 -0.00571429  0.23344785  0.79054264
## Ecatepec de Morelos            0.04146713 -0.00571429  0.08160370  0.16516423
## Ecatzingo                      0.44738058 -0.00571429  0.23344785  0.93776541
## El Oro                        -0.09683642 -0.00571429  0.23344785 -0.18859446
## Emiliano Zapata                0.34151463 -0.00571429  0.15391044  0.88507861
## Gustavo A. Madero              1.89986950 -0.00571429  0.11414173  5.64034689
## Huehuetoca                     0.06620238 -0.00571429  0.23344785  0.14884512
## Hueypoxtla                     0.08465995 -0.00571429  0.31298526  0.16154098
## Huitzilac                     -0.64879100 -0.00571429  0.15391044 -1.63918793
## Huixquilucan                   0.10760397 -0.00571429  0.18572540  0.26294459
## Isidro Fabela                  0.21351661 -0.00571429  0.18572540  0.50870513
## Ixtapaluca                    -0.00068109 -0.00571429  0.15391044  0.01282951
## Ixtapan de la Sal              0.34631182 -0.00571429  0.18572540  0.81684417
## Ixtapan del Oro                0.42327324 -0.00571429  0.31298526  0.76680110
## Ixtlahuaca                     0.24341734 -0.00571429  0.18572540  0.57808702
## Iztacalco                      3.31053141 -0.00571429  0.18572540  7.69504290
## Iztapalapa                     1.98037802 -0.00571429  0.11414173  5.87864445
## Jaltenco                       0.00755830 -0.00571429  0.15391044  0.03383153
## Jantetelco                     0.42730045 -0.00571429  0.23344785  0.89620578
## Jilotepec                      0.31007892 -0.00571429  0.18572540  0.73276907
## Jilotzingo                     0.01164928 -0.00571429  0.18572540  0.04029055
## Jiquipilco                     0.36485488 -0.00571429  0.15391044  0.94457237
## Jiutepec                       0.28604550 -0.00571429  0.23344785  0.60385198
## Jocotitlán                    0.08604397 -0.00571429  0.13118547  0.25333913
## Jojutla                        0.02376588 -0.00571429  0.23344785  0.06101477
## Jonacatepec de Leandro Valle   0.23448221 -0.00571429  0.23344785  0.49713202
## Joquicingo                     0.09431862 -0.00571429  0.15391044  0.25498160
## Juchitepec                    -0.43747374 -0.00571429  0.10088550 -1.35933808
## La Magdalena Contreras         3.80785505 -0.00571429  0.31298526  6.81662978
## La Paz                        -0.02838077 -0.00571429  0.15391044 -0.05777635
## Lerma                          0.01477010 -0.00571429  0.09028051  0.06817512
## Luvianos                       0.10721642 -0.00571429  0.31298526  0.20185992
## Malinalco                      0.27576759 -0.00571429  0.15391044  0.71749089
## Mazatepec                      0.35469114 -0.00571429  0.15391044  0.91866524
## Melchor Ocampo                 0.07343547 -0.00571429  0.31298526  0.14147758
## Metepec                        0.18947611 -0.00571429  0.13118547  0.53890919
## Mexicaltzingo                  0.44736022 -0.00571429  0.31298526  0.80985578
## Miacatlán                     0.32815993 -0.00571429  0.13118547  0.92180705
## Miguel Hidalgo                 2.88083576 -0.00571429  0.13118547  7.96959463
## Milpa Alta                     4.54149654 -0.00571429  0.11414173 13.45931183
## Morelos                        0.33742780 -0.00571429  0.15391044  0.87466137
## Naucalpan de Juárez           0.30194898 -0.00571429  0.11414173  0.91065402
## Nextlalpan                    -0.03645473 -0.00571429  0.11414173 -0.09098878
## Nezahualcóyotl                0.50115932 -0.00571429  0.11414173  1.50029770
## Nicolás Romero                0.12350691 -0.00571429  0.13118547  0.35677210
## Nopaltepec                     0.43630779 -0.00571429  0.94928452  0.45367595
## Ocoyoacac                     -0.00101117 -0.00571429  0.11414173  0.01392076
## Ocuilan                        0.16230965 -0.00571429  0.11414173  0.49733488
## Ocuituco                       0.45881920 -0.00571429  0.23344785  0.96143979
## Otumba                         0.00078201 -0.00571429  0.23344785  0.01344531
## Otzoloapan                     0.30792007 -0.00571429  0.18572540  0.72775965
## Otzolotepec                    0.01036922 -0.00571429  0.15391044  0.04099650
## Ozumba                         0.11314263 -0.00571429  0.13118547  0.32815694
## Papalotla                      0.09712461 -0.00571429  0.31298526  0.18382114
## Polotitlán                    0.59015434 -0.00571429  0.47206007  0.86726512
## Puente de Ixtla               -0.12795927 -0.00571429  0.15391044 -0.31159969
## Rayón                         0.00363773 -0.00571429  0.23344785  0.01935576
## San Antonio la Isla            0.09152224 -0.00571429  0.15391044  0.24785370
## San Felipe del Progreso        0.25841641 -0.00571429  0.15391044  0.67326313
## San José del Rincón          0.45265169 -0.00571429  0.23344785  0.94867496
## San Martín de las Pirámides  0.00009449 -0.00571429  0.18572540  0.01347873
## San Mateo Atenco               0.09725230 -0.00571429  0.31298526  0.18404938
## San Simón de Guerrero         0.02798668 -0.00571429  0.31298526  0.06023936
## Santo Tomás                   0.25920407 -0.00571429  0.31298526  0.47353284
## Soyaniquilpan de Juárez       0.22277735 -0.00571429  0.94928452  0.23451580
## Sultepec                       0.51615315 -0.00571429  0.18572540  1.21094536
## Tecámac                       0.00022609 -0.00571429  0.11414173  0.01758296
## Tejupilco                     -0.06464077 -0.00571429  0.13118547 -0.16269256
## Temamatla                      0.20128425 -0.00571429  0.18572540  0.48032104
## Temascalapa                   -0.03381126 -0.00571429  0.23344785 -0.05815200
## Temascalcingo                  0.05679423 -0.00571429  0.23344785  0.12937319
## Temascaltepec                  0.23541759 -0.00571429  0.10088550  0.75917213
## Temixco                        0.40424972 -0.00571429  0.18572540  0.95128373
## Temoac                         0.63410025 -0.00571429  0.23344785  1.32421702
## Temoaya                        0.18005254 -0.00571429  0.13118547  0.51289127
## Tenancingo                    -0.01925520 -0.00571429  0.15391044 -0.03451549
## Tenango del Aire              -0.05415416 -0.00571429  0.18572540 -0.11240026
## Tenango del Valle              0.00271679 -0.00571429  0.10088550  0.02654414
## Teoloyucan                    -0.00619852 -0.00571429  0.23344785 -0.00100221
## Teotihuacán                  -0.09974196 -0.00571429  0.15391044 -0.23967441
## Tepalcingo                     0.14441834 -0.00571429  0.23344785  0.31072782
## Tepetlaoxtoc                   0.03383955 -0.00571429  0.13118547  0.10920582
## Tepetlixpa                     0.27109306 -0.00571429  0.18572540  0.64230596
## Tepotzotlán                   0.03150028 -0.00571429  0.13118547  0.10274724
## Tepoztlán                    -0.32562726 -0.00571429  0.13118547 -0.88326088
## Tequixquiac                    0.12717342 -0.00571429  0.23344785  0.27503620
## Tetecala                       0.40067870 -0.00571429  0.31298526  0.72641410
## Tetela del Volcán             0.47880979 -0.00571429  0.23344785  1.00281408
## Texcaltitlán                  0.23391857 -0.00571429  0.15391044  0.61081869
## Texcalyacac                   -0.00817712 -0.00571429  0.15391044 -0.00627771
## Texcoco                        0.03057381 -0.00571429  0.09028051  0.12077226
## Tezoyuca                       0.02884594 -0.00571429  0.31298526  0.06177527
## Tianguistenco                 -0.13626139 -0.00571429  0.06301054 -0.52006860
## Timilpan                       0.31583391 -0.00571429  0.15391044  0.81961905
## Tláhuac                       6.24857236 -0.00571429  0.18572540 14.51249650
## Tlalmanalco                    0.10702064 -0.00571429  0.13118547  0.31125448
## Tlalnepantla                  -0.24518797 -0.00571429  0.18572540 -0.55567664
## Tlalnepantla de Baz            0.25768033 -0.00571429  0.13118547  0.72721701
## Tlalpan                        3.53735045 -0.00571429  0.10088550 11.15487509
## Tlaltizapán de Zapata         0.21236172 -0.00571429  0.13118547  0.60209502
## Tlaquiltenango                 0.05996626 -0.00571429  0.15391044  0.16741823
## Tlatlaya                       1.09905005 -0.00571429  0.94928452  1.13389136
## Tlayacapan                    -0.15746976 -0.00571429  0.18572540 -0.35213461
## Toluca                        -0.03825691 -0.00571429  0.10088550 -0.10245618
## Tonanitla                      0.01782151 -0.00571429  0.18572540  0.05461266
## Tonatico                       0.17498718 -0.00571429  0.47206007  0.26300441
## Totolapan                      0.27301347 -0.00571429  0.18572540  0.64676210
## Tultepec                      -0.00134689 -0.00571429  0.18572540  0.01013414
## Tultitlán                    -0.01468854 -0.00571429  0.09028051 -0.02986768
## Valle de Bravo                -0.02563639 -0.00571429  0.13118547 -0.05500377
## Valle de Chalco Solidaridad   -0.53484420 -0.00571429  0.18572540 -1.22779726
## Venustiano Carranza            3.77273368 -0.00571429  0.23344785  7.82021168
## Villa de Allende               0.52745224 -0.00571429  0.23344785  1.10348881
## Villa del Carbón              0.43882908 -0.00571429  0.18572540  1.03152197
## Villa Guerrero                 0.28384665 -0.00571429  0.18572540  0.67189951
## Villa Victoria                 0.43743371 -0.00571429  0.18572540  1.02828412
## Xalatlaco                     -0.03817267 -0.00571429  0.23344785 -0.06717876
## Xochimilco                     9.32722299 -0.00571429  0.18572540 21.65622193
## Xochitepec                     0.40194830 -0.00571429  0.18572540  0.94594350
## Xonacatlán                   -0.00951967 -0.00571429  0.23344785 -0.00787596
## Xoxocotla                      0.23390315 -0.00571429  0.15391044  0.61077938
## Yautepec                       0.26680939 -0.00571429  0.11414173  0.80664418
## Yecapixtla                     0.30399490 -0.00571429  0.10088550  0.97507879
## Zacatepec                     -0.18739238 -0.00571429  0.23344785 -0.37601712
## Zacazonapan                    0.22850884 -0.00571429  0.18572540  0.54349321
## Zacualpan                      0.65163149 -0.00571429  0.23344785  1.36050123
## Zacualpan de Amilpas           0.59318036 -0.00571429  0.23344785  1.23952559
## Zinacantepec                   0.00901322 -0.00571429  0.15391044  0.03754009
## Zumpahuacán                   0.19379577 -0.00571429  0.15391044  0.50854659
## Zumpango                      -0.01007570 -0.00571429  0.10088550 -0.01373133
##                               Pr.z...0.
## Ã\201lvaro Obregón                 0.0000
## Acambay de Ruíz Castañeda      0.2884
## Acolman                          0.4780
## Aculco                           0.1182
## Almoloya de Alquisiras           0.2607
## Almoloya de Juárez              0.3480
## Almoloya del Río                0.0795
## Amacuzac                         0.2489
## Amanalco                         0.5101
## Amatepec                         0.1185
## Amecameca                        0.5376
## Apaxco                           0.3829
## Atenco                           0.4378
## Atizapán                        0.0911
## Atizapán de Zaragoza            0.4698
## Atlacomulco                      0.6003
## Atlatlahucan                     0.2715
## Atlautla                         0.2801
## Axapusco                         0.4150
## Axochiapan                       0.4562
## Ayala                            0.2738
## Ayapango                         0.3694
## Azcapotzalco                     0.0000
## Benito Juárez                   0.0000
## Calimaya                         0.4623
## Capulhuac                        0.4981
## Chalco                           0.0017
## Chapa de Mota                    0.1652
## Chapultepec                      0.3057
## Chiautla                         0.3603
## Chicoloapan                      0.4515
## Chiconcuac                       0.4265
## Chimalhuacán                    0.4876
## Coacalco de Berriozábal         0.4557
## Coatepec Harinas                 0.0879
## Coatetelco                       0.0775
## Coatlán del Río                0.2165
## Cocotitlán                      0.2177
## Coyoacán                        0.0000
## Coyotepec                        0.5328
## Cuajimalpa de Morelos            0.0000
## Cuauhtémoc                      0.0000
## Cuautitlán                      0.6500
## Cuautitlán Izcalli              0.5195
## Cuautla                          0.5512
## Cuernavaca                       0.3736
## Donato Guerra                    0.2146
## Ecatepec de Morelos              0.4344
## Ecatzingo                        0.1742
## El Oro                           0.5748
## Emiliano Zapata                  0.1881
## Gustavo A. Madero                0.0000
## Huehuetoca                       0.4408
## Hueypoxtla                       0.4358
## Huitzilac                        0.9494
## Huixquilucan                     0.3963
## Isidro Fabela                    0.3055
## Ixtapaluca                       0.4949
## Ixtapan de la Sal                0.2070
## Ixtapan del Oro                  0.2216
## Ixtlahuaca                       0.2816
## Iztacalco                        0.0000
## Iztapalapa                       0.0000
## Jaltenco                         0.4865
## Jantetelco                       0.1851
## Jilotepec                        0.2318
## Jilotzingo                       0.4839
## Jiquipilco                       0.1724
## Jiutepec                         0.2730
## Jocotitlán                      0.4000
## Jojutla                          0.4757
## Jonacatepec de Leandro Valle     0.3095
## Joquicingo                       0.3994
## Juchitepec                       0.9130
## La Magdalena Contreras           0.0000
## La Paz                           0.5230
## Lerma                            0.4728
## Luvianos                         0.4200
## Malinalco                        0.2365
## Mazatepec                        0.1791
## Melchor Ocampo                   0.4437
## Metepec                          0.2950
## Mexicaltzingo                    0.2090
## Miacatlán                       0.1783
## Miguel Hidalgo                   0.0000
## Milpa Alta                       0.0000
## Morelos                          0.1909
## Naucalpan de Juárez             0.1812
## Nextlalpan                       0.5362
## Nezahualcóyotl                  0.0668
## Nicolás Romero                  0.3606
## Nopaltepec                       0.3250
## Ocoyoacac                        0.4944
## Ocuilan                          0.3095
## Ocuituco                         0.1682
## Otumba                           0.4946
## Otzoloapan                       0.2334
## Otzolotepec                      0.4836
## Ozumba                           0.3714
## Papalotla                        0.4271
## Polotitlán                      0.1929
## Puente de Ixtla                  0.6223
## Rayón                           0.4923
## San Antonio la Isla              0.4021
## San Felipe del Progreso          0.2504
## San José del Rincón            0.1714
## San Martín de las Pirámides    0.4946
## San Mateo Atenco                 0.4270
## San Simón de Guerrero           0.4760
## Santo Tomás                     0.3179
## Soyaniquilpan de Juárez         0.4073
## Sultepec                         0.1130
## Tecámac                         0.4930
## Tejupilco                        0.5646
## Temamatla                        0.3155
## Temascalapa                      0.5232
## Temascalcingo                    0.4485
## Temascaltepec                    0.2239
## Temixco                          0.1707
## Temoac                           0.0927
## Temoaya                          0.3040
## Tenancingo                       0.5138
## Tenango del Aire                 0.5447
## Tenango del Valle                0.4894
## Teoloyucan                       0.5004
## Teotihuacán                     0.5947
## Tepalcingo                       0.3780
## Tepetlaoxtoc                     0.4565
## Tepetlixpa                       0.2603
## Tepotzotlán                     0.4591
## Tepoztlán                       0.8115
## Tequixquiac                      0.3916
## Tetecala                         0.2338
## Tetela del Volcán               0.1580
## Texcaltitlán                    0.2707
## Texcalyacac                      0.5025
## Texcoco                          0.4519
## Tezoyuca                         0.4754
## Tianguistenco                    0.6985
## Timilpan                         0.2062
## Tláhuac                         0.0000
## Tlalmanalco                      0.3778
## Tlalnepantla                     0.7108
## Tlalnepantla de Baz              0.2335
## Tlalpan                          0.0000
## Tlaltizapán de Zapata           0.2736
## Tlaquiltenango                   0.4335
## Tlatlaya                         0.1284
## Tlayacapan                       0.6376
## Toluca                           0.5408
## Tonanitla                        0.4782
## Tonatico                         0.3963
## Totolapan                        0.2589
## Tultepec                         0.4960
## Tultitlán                       0.5119
## Valle de Bravo                   0.5219
## Valle de Chalco Solidaridad      0.8902
## Venustiano Carranza              0.0000
## Villa de Allende                 0.1349
## Villa del Carbón                0.1511
## Villa Guerrero                   0.2508
## Villa Victoria                   0.1519
## Xalatlaco                        0.5268
## Xochimilco                       0.0000
## Xochitepec                       0.1721
## Xonacatlán                      0.5031
## Xoxocotla                        0.2707
## Yautepec                         0.2099
## Yecapixtla                       0.1648
## Zacatepec                        0.6465
## Zacazonapan                      0.2934
## Zacualpan                        0.0868
## Zacualpan de Amilpas             0.1076
## Zinacantepec                     0.4850
## Zumpahuacán                     0.3055
## Zumpango                         0.5055

7.3.1 Mapping the local Moran’s I

Before mapping the local Moran’s I map, it is wise to append the local Moran’s I dataframe (i.e. localMI) onto hunan SpatialPolygonDataFrame. The code chunks below can be used to perform the task. The out SpatialPolygonDataFrame is called mexi_combined_filtered.localMI.

mexi_combined_filtered.localMI <- cbind(mexi_combined_filtered,localMI)

7.3.2 Mapping local Moran’s I values

Using choropleth mapping functions of tmap package, to plot the local Moran’s I values by using the code chinks below.

tm_shape(mexi_combined_filtered.localMI) +
  tm_fill(col = "Ii", 
          style = "pretty",
          palette = "RdBu",
          title = "local moran statistics") +
  tm_borders(alpha = 0.5)

### Mapping local Moran’s I p-values The choropleth shows there is evidence for both positive and negative Ii values. However, it is useful to consider the p-values for each of these values, as consider above.

The code chunks below produce a choropleth map of Moran’s I p-values by using functions of tmap package.

tm_shape(mexi_combined_filtered.localMI) +
  tm_fill(col = "Pr.z...0.", 
          breaks=c(-Inf, 0.001, 0.01, 0.05, 0.1, Inf),
          palette="-Blues", 
          title = "local Moran's I p-values") +
  tm_borders(alpha = 0.5)

Mapping both local Moran’s I values and p-values For effective interpretation, it is better to plot both the local Moran’s I values map and its correpence p-values map next to each other.

The code chunk below will be used to create such visualisation.

localMI.map <- tm_shape(mexi_combined_filtered.localMI) +
  tm_fill(col = "Ii", 
          style = "pretty", 
          title = "local moran statistics") +
  tm_borders(alpha = 0.5)

pvalue.map <- tm_shape(mexi_combined_filtered.localMI) +
  tm_fill(col = "Pr.z...0.", 
          breaks=c(-Inf, 0.001, 0.01, 0.05, 0.1, Inf),
          palette="-Blues", 
          title = "local Moran's I p-values") +
  tm_borders(alpha = 0.5)

tmap_arrange(localMI.map, pvalue.map, asp=1, ncol=2)

## LISA Cluster Map

The LISA Cluster Map shows the significant locations color coded by type of spatial autocorrelation. The first step before we can generate the LISA cluster map is to plot the Moran scatterplot.

Plotting Moran scatterplot The Moran scatterplot is an illustration of the relationship between the values of the chosen attribute at each location and the average value of the same attribute at neighboring locations.

nci <- moran.plot(mexi_combined_filtered$rate_of_covid, rswm_q_mexi_combined_filtered, labels=as.character(mexi_combined_filtered$NOMGEO), xlab="COVID 19 Rate (Pert 10,000)", ylab="Spatially Lag COVID 19 Rate")

The pattern could be observed where the higher COVID 19 rate is being clustered together with the high COVID 19 rate area, while the lower rate are directly in the opposite.

7.3.3 Plotting Moran scatterplot with standardised variable

Scaling the variable of interest First we will use scale() to centers and scales the variable. here centering is done by subtracting the mean (omitting NAs) the corresponding columns, and scaling is done by dividing the (centered) variable by their standard deviations.

mexi_combined_filtered$Z.rate_of_covid <- scale(mexi_combined_filtered$rate_of_covid) %>% as.vector 

The as.vector() added to the end is to make sure that the data type we get out of this is a vector, that map neatly into the dataframe.

Now, we are ready to plot the Moran scatterplot again by using the code chunk below.

nci2 <- moran.plot(mexi_combined_filtered$Z.rate_of_covid, rswm_q_mexi_combined_filtered, labels=as.character(mexi_combined_filtered$NOMGEO), xlab="z-rate_of_COVID19 (Per 10,000)", ylab="Spatially Lag z-rate_of_COVID19")

The code chunks below show the steps to prepare a LISA cluster map. Next, we centers the variable of interest around its mean. This is follow by centering the local Moran’s around the mean. Next, we will set a statistical significance level for the local Moran. These four command lines define the high-high, low-low, low-high and high-low categories. Lastly, places non-significant Moran in the category 0.

quadrant <- vector(mode="numeric",length=nrow(localMI))
DV <- mexi_combined_filtered$rate_of_covid - mean(mexi_combined_filtered$rate_of_covid)
C_mI <- localMI[,1] - mean(localMI[,1])    
signif <- 0.05       

quadrant[DV >0 & C_mI>0] <- 4      
quadrant[DV <0 & C_mI<0] <- 1      
quadrant[DV <0 & C_mI>0] <- 2
quadrant[DV >0 & C_mI<0] <- 3
quadrant[localMI[,5]>signif] <- 0

Now, we can build the LISA map by using the code chunks below.

mexi_combined_filtered.localMI$quadrant <- quadrant
colors <- c("#ffffff", "#2c7bb6", "#abd9e9", "#fdae61", "#d7191c")
clusters <- c("insignificant", "low-low", "low-high", "high-low", "high-high")

tm_shape(mexi_combined_filtered.localMI) +
  tm_fill(col = "quadrant", 
          style = "cat", 
          palette = colors[c(sort(unique(quadrant)))+1], 
          labels = clusters[c(sort(unique(quadrant)))+1]) +
  tm_view(set.zoom.limits = c(11,17)) +
  tm_borders(alpha=0.5)

# Hot Spot and Cold Spot Area Analysis

Beside detecting cluster and outliers, localised spatial statistics can be also used to detect hot spot and/or cold spot areas.

The term ‘hot spot’ has been used generically across disciplines to describe a region or value that is higher relative to its surroundings (Lepers et al 2005, Aben et al 2012, Isobe et al 2015).

7.4 Getis and Ord’s G-Statistics

An alternative spatial statistics to detect spatial anomalies is the Getis and Ord’s G-statistics (Getis and Ord, 1972; Ord and Getis, 1995). It looks at neighbours within a defined proximity to identify where either high or low values clutser 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.

The analysis consists of three steps:

  • Deriving spatial weight matrix
  • Computing Gi statistics
  • Mapping Gi statistics

7.5 Deriving distance-based weight matrix

First, we need to define a new set of neighbours. Whist the spatial autocorrelation considered units which shared borders, for Getis-Ord we are defining neighbours based on distance or also known as proximity.

There are two type of distance-based proximity matrix, they are:

  • fixed distance weight matrix; and
  • adaptive distance weight matrix.

7.6 Creating adaptive proximity matrix

But since fixed distance weight matrix is not so effective as seen or demonstrate from the Local Moran’s I result. Therefore, we will proceed to use adaptive approach.

Instead of using a fixed distance to deirved the proximity matrix, a fixed number of neighbours can be also used to derived the proximity matrix. This is called adaptive distance method.

The code chunk below will be used to derive an adaptive distance weight matrix.

coords <- coordinates(mexi_combined_filtered)
#knb <- knn2nb(knearneigh(coords, k=8, longlat = TRUE), row.names=row.names(mexi_combined_filtered$rate_of_covid))
knn6_lw <- nb2listw(knn6, style = 'B')
summary(knn6_lw)
## Characteristics of weights list object:
## Neighbour list object:
## Number of regions: 176 
## Number of nonzero links: 1056 
## Percentage nonzero weights: 3.409091 
## Average number of links: 6 
## Non-symmetric neighbours list
## Link number distribution:
## 
##   6 
## 176 
## 176 least connected regions:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 with 6 links
## 176 most connected regions:
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 with 6 links
## 
## Weights style: B 
## Weights constants summary:
##     n    nn   S0   S1    S2
## B 176 30976 1056 1874 26004

Next, we will plot the centroids of the counties and neighbours by using the code chunk below.

plot(mexi_combined_filtered, border="lightgrey")
plot(knn6, coords, pch = 19, cex = 0.6, add = TRUE, col = "red")

In the choropleth map above, counties shaded in red are the hot spot areas and counties shaded in blue are the cold spot areas. The darkness of the colours representing the intensity of the Gi values.

Useful to conduct this since the Local Moran’s I does not help us identify the degree of clustering that is present. We use the adaptive distance weights matrix as the fixed distance one we created is not meaningful.

Gi statistics using adaptive distance

The code chunk below are used to compute the Gi values for Rate of COVID 19 by using an adaptive distance weight matrix (i.e knb_lw).

fips <- order(mexi_combined_filtered$NOMGEO)
gi.adaptive <- localG(mexi_combined_filtered$rate_of_covid, knn6_lw)
mexi_combined_filtered.gi <- cbind(mexi_combined_filtered, as.matrix(gi.adaptive))
names(mexi_combined_filtered.gi)[95] <- "gstat_adaptive"
tm_shape(mexi_combined_filtered.gi) +
  tm_fill(col = "gstat_adaptive", 
          style = "pretty",
          palette="-RdBu",
          title = "local Gi") +
  tm_borders(alpha = 0.5) + tm_text("NOMGEO", size=0.2) 

As we can observed the result above, the hotspot are being highlighted in red and the cold spots have a color marked as blue. The color describe the density and the weight of the Gi values. Which allows us to better analyse the Moran’s I method. The area of Mexico States and City are being highlighted with red as it has a considerably high clustered population than the Morelos State area. But looking at Mexico States, the area are wide apart with the COVID spread around the corner while the critical zone were near the Mexico City. And even though Mexico City is small but the population are higher due to the area in measurement is smaller compared to Mexico States. Therefor it greatly affect the influence of spreading the COVID.

8 Conclusion

At last, we are able to identify the result from the outcome after analyze the Study Area with various techniques and methods like EDA, Thematic mapping, Weight Matrix and spatial analysis. With a higher intense COVID 19 range of 109.3 to 170.4 in Mexico City and have lesser municipalities as compared to Mexico State and Morelos State. It could be due to higher densely population and are more closely connected that increase the population to falls into COVID 19.