Создайте карту пространственных кластеров LISA в R

Question

Я хотел бы создать карту, показывающую локальный пространственный кластер явления, предпочтительно используя Local Moran (LISA).

В примере, воспроизводимом ниже, я вычисляю индекс локального морана, используя spdep, но я хотел бы знать, есть ли простой способ сопоставить кланы, предварительно используя ggplot2. Помогите ?

library(UScensus2000tract) 
library(ggplot2) 
library(spdep) 

# load data 
data("oregon.tract") 

# plot Census Tract map 
plot(oregon.tract) 

# create Queens contiguity matrix 
spatmatrix <- poly2nb(oregon.tract) 

#calculate the local moran of the distribution of black population 
lmoran <- localmoran(oregon.tract@data$black, nb2listw(spatmatrix)) 

Теперь, чтобы сделать этот пример более похож на мои реальные данные, у меня есть некоторые NA значения в моем файле формы, которые представляют собой отверстие в многоугольнике, поэтому эти участки не должны использоваться в расчетах.

oregon.tract@data$black[3:5] <- NA 

Answer 1

Я не думаю, что этот ответ достоин щедрости, но, возможно, он приблизит вас к ответу. Поскольку я ничего не знаю о localmoran, я только догадывался о заполнении.

library(UScensus2000tract) 
library(ggplot2) 
library(spdep) 

# load data 
data("oregon.tract") 

# plot Census Tract map 
plot(oregon.tract) 

# create Queens contiguity matrix 
spatmatrix <- poly2nb(oregon.tract) 

#calculate the local moran of the distribution of black population 
lmoran <- localmoran(oregon.tract@data$black, nb2listw(spatmatrix)) 

# get our id from the rownames in a data.frame 
oregon.tract@data$id <- rownames(oregon.tract@data) 
oregon.tract@data$lmoran_ii <- lmoran[,1] 
oregon_df <- merge(
    # convert to a data.frame 
    fortify(oregon.tract, region="id"), 
    oregon.tract@data, 
    by="id" 
) 

ggplot(data=oregon_df, aes(x=long,y=lat,group=group)) + 
    geom_polygon(fill=scales::col_numeric("Blues",domain=c(-1,5))(oregon_df$lmoran_ii)) + 
    geom_path(color="white") 

Answer 2

Вот стратегия:

library(UScensus2000tract) 
library(spdep) 
library(ggplot2) 
library(dplyr) 

# load data 
data("oregon.tract") 
# plot Census Tract map 
plot(oregon.tract) 

# create Queens contiguity matrix 
spatmatrix <- poly2nb(oregon.tract) 

# create a neighbours list with spatial weights 
listw <- nb2listw(spatmatrix) 

# calculate the local moran of the distribution of white population 
lmoran <- localmoran(oregon.tract$white, listw) 
summary(lmoran) 

# padronize the variable and save it to a new column 
oregon.tract$s_white <- scale(oregon.tract$white) %>% as.vector() 

# create a spatially lagged variable and save it to a new column 
oregon.tract$lag_s_white <- lag.listw(listw, oregon.tract$s_white) 

# summary of variables, to inform the analysis 
summary(oregon.tract$s_white) 
summary(oregon.tract$lag_s_white) 

# moran sccaterplot, in basic graphics (with identification of influential observations) 
x <- oregon.tract$s_white 
y <- oregon.tract$lag_s_white %>% as.vector() 
xx <- data_frame(x, y) 

moran.plot(x, listw) 

# moran sccaterplot, in ggplot 
# (without identification of influential observations - which is possible but requires more effort) 
ggplot(xx, aes(x, y)) + geom_point() + geom_smooth(method = 'lm', se = F) + geom_hline(yintercept = 0, linetype = 'dashed') + geom_vline(xintercept = 0, linetype = 'dashed') 

# create a new variable identifying the moran plot quadrant for each observation, dismissing the non-significant ones 
oregon.tract$quad_sig <- NA 

# high-high quadrant 
oregon.tract[(oregon.tract$s_white >= 0 & 
       oregon.tract$lag_s_white >= 0) & 
       (lmoran[, 5] <= 0.05), "quad_sig"] <- "high-high" 
# low-low quadrant 
oregon.tract[(oregon.tract$s_white <= 0 & 
       oregon.tract$lag_s_white <= 0) & 
       (lmoran[, 5] <= 0.05), "quad_sig"] <- "low-low" 
# high-low quadrant 
oregon.tract[(oregon.tract$s_white >= 0 & 
       oregon.tract$lag_s_white <= 0) & 
       (lmoran[, 5] <= 0.05), "quad_sig"] <- "high-low" 
# low-high quadrant 
oregon.tract@data[(oregon.tract$s_white <= 0 
       & oregon.tract$lag_s_white >= 0) & 
       (lmoran[, 5] <= 0.05), "quad_sig"] <- "low-high" 
# non-significant observations 
oregon.tract@data[(lmoran[, 5] > 0.05), "quad_sig"] <- "not signif." 

oregon.tract$quad_sig <- as.factor(oregon.tract$quad_sig) 
oregon.tract@data$id <- rownames(oregon.tract@data) 

# plotting the map 
df <- fortify(oregon.tract, region="id") 
df <- left_join(df, oregon.tract@data) 
df %>% 
    ggplot(aes(long, lat, group = group, fill = quad_sig)) + 
    geom_polygon(color = "white", size = .05) + coord_equal() + 
    theme_void() + scale_fill_brewer(palette = "Set1") 

Этот ответ был основан на this page, предложенный Eli Knaap on twitter, а также заимствованные из ответа по @timelyportfolio на этот вопрос.

Я использовал переменную white вместо black, потому что black имел менее явные результаты.

Относительно Nas, localmoran() включает аргумент na.action, о котором говорится в документации:

na.action является функцией (na.fail по умолчанию), а также может быть na.omit или> na.exclude - в этих случаях список весов будет подмножеством для удаления NA в данных. Возможно, необходимо установить для параметра zero.policy значение TRUE, потому что это подмножество может создавать наблюдения соседей. Обратите внимание, что только списки весов, созданные без использования аргумента glist для nb2listw, могут быть подмножествами. Если используется na.pass, то при вычислении пространственного запаздывания вместо значений NA заменяется ноль.

Я пробовал:

oregon.tract@data$white[3:5] <- NA 
lmoran <- localmoran(oregon.tract@data$white, listw, zero.policy = TRUE, 
       na.action = na.exclude) 

Но столкнулись с проблемами в lag.listw и не было времени, чтобы посмотреть на него. Сожалею.