在 R 中创建空间集群 LISA 的地图
Create a map of spatial clusters LISA in R
我想创建一个显示现象的局部空间聚类的地图，最好使用 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
我认为这个答案不值得悬赏，但也许它会让您更接近答案。由于我对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")
这是一个策略：
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 scatterplot, 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, suggested by Eli Knaap on twitter，也借鉴了@timelyportfolio对该问题的回答。
我使用变量 white 而不是 black 因为 black 的结果不太明确。
关于 NAs，localmoran() 包括参数 na.action，关于它的文档说：
na.action is a function (default na.fail), can also be na.omit or > na.exclude - in these cases the weights list will be subsetted to remove NAs in the data. It may be necessary to set zero.policy to TRUE because this subsetting may create no-neighbour observations. Note that only weights lists created without using the glist argument to nb2listw may be subsetted. If na.pass is used, zero is substituted for NA values in calculating the spatial lag.
我试过了：
oregon.tract@data$white[3:5] <- NA
lmoran <- localmoran(oregon.tract@data$white, listw, zero.policy = TRUE,
na.action = na.exclude)
但是 运行 遇到了 lag.listw 中的问题，但没有时间去研究它。对不起。
这显然已经很晚了，但我在做类似的事情时遇到了 post。这使用了 rgeoda package，当问题被 posted 时它还没有出现，但它是由 GeoDa 人员开发的，用于将该软件的一些功能移植到 R。sf 有与此同时，它也真正起飞了，这使得操作空间数据变得非常容易； rgeoda 函数通常需要 sf 个对象。
与另一个 poster 一样，我使用的是白人而不是黑人，因为集群显示得更好。我将原始数据（缺少少数观察值）转换为 sf。 rgeoda::local_moran 在缺少数据时不起作用，但如果您复制并删除了缺少的观察值，则可以 运行 分析并通过 ID 将它们重新组合在一起。使用右连接，以便保留原始数据中的所有 ID，包括缺失值。
因为它模仿 GeoDa，所以 LISA 对象中存储了与 local_moran returns 相同的颜色和标签。提取它们并将它们用作调色板。由于调色板已命名，并且这些名称不包含“NA”，因此您可以将 NA 值添加到调色板矢量，或手动为 NA 值指定颜色以确保仍绘制这些形状。我把它变成了绿色，这样它就可以看到了（左上角）。
library(UScensus2000tract)
library(ggplot2)
library(dplyr)
library(sf)
library(rgeoda)
# load data
data("oregon.tract")
oregon.tract@data$white[3:5] <- NA
ore_sf <- st_as_sf(oregon.tract) %>%
tibble::rownames_to_column("id")
to_clust <- ore_sf %>%
filter(!is.na(white))
queen_wts <- queen_weights(to_clust)
moran <- local_moran(queen_wts, st_drop_geometry(to_clust["white"]))
moran_lbls <- lisa_labels(moran)
moran_colors <- setNames(lisa_colors(moran), moran_lbls)
ore_clustered <- to_clust %>%
st_drop_geometry() %>%
select(id) %>%
mutate(cluster_num = lisa_clusters(moran) + 1, # add 1 bc clusters are zero-indexed
cluster = factor(moran_lbls[cluster_num], levels = moran_lbls)) %>%
right_join(ore_sf, by = "id") %>%
st_as_sf()
ggplot(ore_clustered, aes(fill = cluster)) +
geom_sf(color = "white", size = 0) +
scale_fill_manual(values = moran_colors, na.value = "green") +
theme_dark()
我想创建一个显示现象的局部空间聚类的地图，最好使用 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
我认为这个答案不值得悬赏，但也许它会让您更接近答案。由于我对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")
这是一个策略：
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 scatterplot, 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, suggested by Eli Knaap on twitter，也借鉴了@timelyportfolio对该问题的回答。
我使用变量 white 而不是 black 因为 black 的结果不太明确。
关于 NAs，localmoran() 包括参数 na.action，关于它的文档说：
na.action is a function (default na.fail), can also be na.omit or > na.exclude - in these cases the weights list will be subsetted to remove NAs in the data. It may be necessary to set zero.policy to TRUE because this subsetting may create no-neighbour observations. Note that only weights lists created without using the glist argument to nb2listw may be subsetted. If na.pass is used, zero is substituted for NA values in calculating the spatial lag.
我试过了：
oregon.tract@data$white[3:5] <- NA
lmoran <- localmoran(oregon.tract@data$white, listw, zero.policy = TRUE,
na.action = na.exclude)
但是 运行 遇到了 lag.listw 中的问题，但没有时间去研究它。对不起。
这显然已经很晚了，但我在做类似的事情时遇到了 post。这使用了 rgeoda package，当问题被 posted 时它还没有出现，但它是由 GeoDa 人员开发的，用于将该软件的一些功能移植到 R。sf 有与此同时，它也真正起飞了，这使得操作空间数据变得非常容易； rgeoda 函数通常需要 sf 个对象。
与另一个 poster 一样，我使用的是白人而不是黑人，因为集群显示得更好。我将原始数据（缺少少数观察值）转换为 sf。 rgeoda::local_moran 在缺少数据时不起作用，但如果您复制并删除了缺少的观察值，则可以 运行 分析并通过 ID 将它们重新组合在一起。使用右连接，以便保留原始数据中的所有 ID，包括缺失值。
因为它模仿 GeoDa，所以 LISA 对象中存储了与 local_moran returns 相同的颜色和标签。提取它们并将它们用作调色板。由于调色板已命名，并且这些名称不包含“NA”，因此您可以将 NA 值添加到调色板矢量，或手动为 NA 值指定颜色以确保仍绘制这些形状。我把它变成了绿色，这样它就可以看到了（左上角）。
library(UScensus2000tract)
library(ggplot2)
library(dplyr)
library(sf)
library(rgeoda)
# load data
data("oregon.tract")
oregon.tract@data$white[3:5] <- NA
ore_sf <- st_as_sf(oregon.tract) %>%
tibble::rownames_to_column("id")
to_clust <- ore_sf %>%
filter(!is.na(white))
queen_wts <- queen_weights(to_clust)
moran <- local_moran(queen_wts, st_drop_geometry(to_clust["white"]))
moran_lbls <- lisa_labels(moran)
moran_colors <- setNames(lisa_colors(moran), moran_lbls)
ore_clustered <- to_clust %>%
st_drop_geometry() %>%
select(id) %>%
mutate(cluster_num = lisa_clusters(moran) + 1, # add 1 bc clusters are zero-indexed
cluster = factor(moran_lbls[cluster_num], levels = moran_lbls)) %>%
right_join(ore_sf, by = "id") %>%
st_as_sf()
ggplot(ore_clustered, aes(fill = cluster)) +
geom_sf(color = "white", size = 0) +
scale_fill_manual(values = moran_colors, na.value = "green") +
theme_dark()