【问题标题】:Neighbour cut off邻居断了
【发布时间】:2023-03-03 22:41:02
【问题描述】:

我正在尝试根据某个国家/地区的省份是否与感兴趣的国家/地区接壤来对空间多边形数据框进行子集化。

假设我对捷克* (CZE) 感兴趣,我的分析单位是第一级行政区划(见图表,蓝色的 CZE)。

我如何对数据进行子集化以包括所有捷克省份以及其他国家的毗邻省份,例如巴伐利亚?

code for data frame and plot

【问题讨论】:

    标签: r maps spatial nearest-neighbor


    【解决方案1】:

    有趣的问题,这是我根据this 使用poly2nb 提出的。我以前没有用过这个,所以它可能有点笨拙,但至少它似乎有效并且可能会让你上路。在this 线程中给出了另一个选项,使用gRelate

    library(spdep)
    
    countries<-c("CZE","DEU","AUT","SVK","POL")
    EUR<-getCountries(countries,level=1) # Level 1
    CZE<-subset(EUR,ISO=="CZE")
    plot(EUR)
    
    # create neighbour structure with ID_1 as names
    # each item in this list contains id of all neighbouring polygons
    nb = poly2nb(pl=EUR, row.names=EUR$ID_1)
    
    # plot nb structure, collection of lines and points
    plot(nb, coordinates(EUR), col=2, lty=2, add=T)
    
    # get ID_1 names from nb items, stored in attributes
    nb_attr = attributes(nb)$region.id
    
    # take only those nb items that belong to the CZE polygons
    nb1 = which(nb_attr %in% CZE$ID_1)
    
    # convert the numeric id in the selected items to ID_1 format
    neighbour_ids = unique(nb_attr[unlist(nb[nb1])])
    
    # select polygons in EUR dataset that correspond with the found ID_1's 
    CZE_neighbours = subset(EUR,ID_1 %in% neighbour_ids)
    
    # plot neighbours and CZE
    plot(CZE_neighbours, col=3, add=T)
    plot(CZE, col='steelblue4', add=T)
    
    # add legend
    legend('topleft',fill=c('steelblue4','green'),legend=c('CZE','neighbours'))
    

    为了完整起见,在rgeos 包中使用gRelate 的替代方法:

    # get relations between EUR and CZE polygons
    x<-gRelate(EUR, CZE, byid = TRUE)
    
    # Select polygon ids which border at least 1 of the CZE polygons 
    # ("FF2F11212" are bordering polygons) 
    id_bordering = EUR$ID_1[apply(x,2,function(y) sum(which(y=="FF2F11212"))>=1)]
    poly_bordering = subset(EUR, ID_1 %in% id_bordering)
    
    # plot results
    plot(EUR)
    plot(poly_bordering, add=T, col=4)
    

    【讨论】: