15.2 Global Measures

Average level of spatial autocorrelation across all observations

joincount.test is similar to Moran’s \(I\) in {spdep}

The BB join count statistic for k-coloured factors count test for spatial autocorrelation using a spatial weights matrix in weights list form for testing whether same-colour joins occur more frequently than would be expected if the zones were labelled in a spatially random way.

?joincount.test
?joincount.multi
args(joincount.test)
## function (fx, listw, zero.policy = attr(listw, "zero.policy"), 
##     alternative = "greater", sampling = "nonfree", spChk = NULL, 
##     adjust.n = TRUE) 
## NULL
args(moran.test)
## function (x, listw, randomisation = TRUE, zero.policy = attr(listw, 
##     "zero.policy"), alternative = "greater", rank = FALSE, na.action = na.fail, 
##     spChk = NULL, adjust.n = TRUE, drop.EI2 = FALSE) 
## NULL
Types <- pol_pres15 |> 
        st_drop_geometry() |> 
        subset(select = types, drop = TRUE)
Types |> 
    table()
## Types
##          Rural          Urban    Urban/rural Warsaw Borough 
##           1563            303            611             18
Types |> joincount.multi(listw = lw_q_B) %>%
  head
##                               Joincount     Expected     Variance     z-value
## Rural:Rural                        3087 2793.9201781 1126.5342033   8.7320000
## Urban:Urban                         110  104.7185351   93.2993687   0.5467831
## Urban/rural:Urban/rural             656  426.5255306  331.7590322  12.5986206
## Warsaw Borough:Warsaw Borough        41    0.3501833    0.3474277  68.9646203
## Urban:Rural                         668 1083.9408630  708.2086432 -15.6297121
## Urban/rural:Rural                  2359 2185.7685388 1267.1313345   4.8664913

Using an inverse distance based listw object releases different results. We first need to identify the Neighbourhood contiguity by distance or the neighbours of region points calculated by Euclidean distance in the metric of the points between the lower and the upper bounds.

?dnearneigh
  • lower distance bound: 0
  • upper distance bound: 18300km
nb_d183 <- coords |> dnearneigh(0, 18300)
nb_d183
## Neighbour list object:
## Number of regions: 2495 
## Number of nonzero links: 21086 
## Percentage nonzero weights: 0.3387296 
## Average number of links: 8.451303

Then consider the Spatial link distance measures from the neighbour list object (nb_d183).

gwts <- nb_d183 |> 
    nbdists(coords) |> 
    lapply(function(x) 1/(x/1000))

Calculate again the Spatial weights for neighbours lists with the nb2listw function.

lw_d183_idw_B <- nb_d183 |> 
  nb2listw(glist=gwts, style="B")
Types |> 
  joincount.multi(listw = lw_d183_idw_B) 
##                                Joincount   Expected   Variance  z-value
## Rural:Rural                   3.4648e+02 3.6123e+02 4.9314e+01  -2.1003
## Urban:Urban                   2.9045e+01 1.3539e+01 2.2281e+00  10.3877
## Urban/rural:Urban/rural       4.6498e+01 5.5145e+01 9.6134e+00  -2.7891
## Warsaw Borough:Warsaw Borough 1.6822e+01 4.5275e-02 6.6083e-03 206.3805
## Urban:Rural                   2.0206e+02 1.4014e+02 2.3645e+01  12.7338
## Urban/rural:Rural             2.2517e+02 2.8260e+02 3.5892e+01  -9.5852
## Urban/rural:Urban             3.6499e+01 5.4784e+01 8.8586e+00  -6.1434
## Warsaw Borough:Rural          5.6502e+00 8.3253e+00 1.7260e+00  -2.0362
## Warsaw Borough:Urban          9.1800e+00 1.6139e+00 2.5392e-01  15.0150
## Warsaw Borough:Urban/rural    3.2676e+00 3.2545e+00 5.5180e-01   0.0177
## Jtot                          4.8183e+02 4.9072e+02 4.1570e+01  -1.3781
I_turnout <- pol_pres15 |> 
        st_drop_geometry() |> 
        subset(select = I_turnout, drop = TRUE)
options(scipen = 999)
I_turnout |> 
  moran.test(listw = lw_q_B, randomisation = FALSE) 
## 
##  Moran I test under normality
## 
## data:  I_turnout  
## weights: lw_q_B    
## 
## Moran I statistic standard deviate = 58.461, p-value <
## 0.00000000000000022
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.6914339743     -0.0004009623      0.0001400449

lm.morantest extracts the residuals used for testing to compare with the standard test.

lm(I_turnout ~ 1, pol_pres15) |> 
    lm.morantest(listw = lw_q_B) 
## 
##  Global Moran I for regression residuals
## 
## data:  
## model: lm(formula = I_turnout ~ 1, data = pol_pres15)
## weights: lw_q_B
## 
## Moran I statistic standard deviate = 58.461, p-value <
## 0.00000000000000022
## alternative hypothesis: greater
## sample estimates:
## Observed Moran I      Expectation         Variance 
##     0.6914339743    -0.0004009623     0.0001400449
mtr <- I_turnout |> 
    moran.test(listw = lw_q_B) 

mtr |> 
    glance_htest()
## Moran I statistic       Expectation          Variance       Std deviate 
##      0.6914339743     -0.0004009623      0.0001400522     58.4598351527 
##           p.value 
##      0.0000000000

In the early 1970s, interest was shown in Monte Carlo tests, also known as Hope-type tests and as permutation bootstrap.

Permutation test for Moran’s I statistic is calculated by using nsim random permutations of x for the given spatial weighting scheme, to establish the rank of the observed statistic in relation to the nsim simulated values.

?moran.mc
set.seed(1)
mmc <- I_turnout |> 
    moran.mc(listw = lw_q_B, 
             nsim = 999, 
             return_boot = TRUE) 
c("Permutation bootstrap" = var(mmc$t), 
  "Analytical randomisation" = unname(mtr$estimate[3]))
##    Permutation bootstrap Analytical randomisation 
##             0.0001441596             0.0001400522