合并两个空间对象的 ID 以创建等值线图

Merging IDs of two spatial objects to create a choropleth map

目标: 根据加泰罗尼亚各市使用人口创建等值线图。

可重现数据: 好的,我的第一步是下载人口和自治市 shapefile。

人口:https://www.idescat.cat/cat/idescat/biblioteca/docs/publicacions/gridpoblacio01012016.zip 自治市边界:http://auriga.icc.cat/bseccen_etrs89/bseccenv10sh1f1_2002a2016_0.zip

目前的步数:

都导入了,给了他们相同的坐标:

catapop<-readOGR("location","rp2016_qtree_level2_ofus_allvar")
catasense<-readOGR("location","bseccenv10sh1f1_20160101_0")

catapop<-spTransform(catapop,CRSobj = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs+towsgs84=0,0,0")
catasense<-spTransform(catasense,CRSobj = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs+towsgs84=0,0,0")

问题

当我查看 shapefile 中的数据时,它们包含我需要的内容:

但是,当我比较 catasense 中的 ID 与 catapop 中的 ID 时,我不知道发生了什么,也不知道如何将 catapop 的 ID 与 catasense 相匹配。

我想保留 catasense "MUNICIPI" 的 ID,因为它们似乎是加泰罗尼亚 public 数据中最标准的。

任何关于如何匹配 ID 并为人口 "TOTAL" 创建叶绿体图的任何想法都将不胜感激!

如果需要任何说明,请告诉我!

所以第一步是将这两个表导入到R

tab1

tab2

library(readxl) 
pop <- read_excel("Downloads/rp2016/pop.xlsx")  
cod <- read_excel("Downloads/rp2016/cod.xlsx")

names(cod) <- c("Codi", "Nom2", "Codi comarca", "Nom comarca") 
codf <- merge(cod, pop, by.x = "Nom2", by.y = "Nom", all.x = TRUE)

#I make a treatment in the Codi field to put 0 in front of a code that starts with 8.

b <- codf$Codi 

b[grep("^8", b)] <- paste0("0",b[grep("^8", b)]) 

codf$Codi <- b


data2 <- catasense@data 

codf2 <- merge(data2, codf, by.x = "MUNICIPI", by.y = "Codi", all.x = TRUE, sort = FALSE)
catasense@data$pop <- codf2$`Població (2016)` 
catasense@data$name <- codf2$Nom2

library(leaflet) 


pal <- colorNumeric("viridis", NULL)


map <- leaflet(catasense) %>% 
  addPolygons(color = "#444444", weight = 1, smoothFactor = 0.5,
              opacity = 1.0, fillOpacity = 1,
              fillColor = ~pal(log10(as.numeric(pop))),
              popup = ~paste0("<b>", name, "</b>", " <br> ", "pop:", pop, "<br>"
              ),
              label = ~paste0(name),
              highlightOptions = highlightOptions(color = "white", weight = 2,
                                              bringToFront = TRUE)) %>%
   addLegend(pal = pal, values = ~log10(pop), opacity = 1.0,
        labFormat = labelFormat(transform = function(x) round(10^x)))


map

#you can save leaflet map in html
library(htmlwidgets) 
saveWidget(map, file="cata2.html")

download and open this html file and see the map