0

R 新手，我有一个很长的问题：

我有一个 shapefile/地图，我的目标是根据该多边形和与其相邻的每个多边形的属性为该地图中的每个多边形计算某个索引。

我有一个邻接矩阵——我认为它与“一阶皇后邻接权重矩阵”相同，尽管我不确定——它描述了哪些多边形与其他多边形接壤，例如，

POLYID A B C D E  
    A  0 0 1 0 1  
    B  0 0 1 0 0    
    C  1 1 0 1 0     
    D  0 0 1 0 1     
    E  1 0 0 1 0

例如，上面表示多边形“C”和“E”与多边形“A”相邻；多边形“B”仅邻接多边形“C”等。

我拥有的属性表每行有一个多边形：

POLYID TOT L10K 10_15K 15_20K ...  
     A 500   24     30     77 ...

其中 TOT、L10K 等是我用来计算指数的变量。

我的数据中有 525 个多边形/行，所以我想使用邻接矩阵来确定将哪些行的属性合并到感兴趣索引的计算中。现在，当我对与相邻多边形的一个“束”相对应的行进行子集化时，我可以计算索引，然后使用循环（如果感兴趣，我正在计算百分位差距指数，一种衡量当地收入隔离的方法）。例如，子集底特律市学校的“社区”：

Detroit <- UNSD00[c(142,150,164,221,226,236,295,327,157,177,178,364,233,373,418,424,449,451,487),]

然后记录边缘列的比例和运行总计：

catprops <- vector()
for(i in 4:19)
{
  catprops[(i-3)]<-sum(Detroit[,i])/sum(Detroit[,3])
}
catprops <- as.data.frame(catprops)
catprops[,2]<-cumsum(catprops[,1])

列 4:19 是属性表中的必要列。

然后我使用下面的代码来计算索引——注意循环有“i in 1:19”，因为底特律子集有 19 个多边形。

cgidistsum <- 0
for(i in 1:19)
{  
   pranks <- vector()
   for(j in 4:19)
    {
      if (Detroit[i,j]==0)
        pranks <- append(pranks,0)
      else if (j == 4)
      pranks <- append(pranks,seq(0,catprops[1,2],by=catprops[1,2]/Detroit[i,j]))
      else 
        pranks <- append(pranks,seq(catprops[j-4,2],catprops[j-3,2],by=catprops[j-3,1]/Detroit[i,j]))
    }
  distpranks <- vector()
  distpranks<-abs(pranks-median(pranks))
  cgidistsum <- cgidistsum + sum(distpranks)
  }
cgi <- (.25-(cgidistsum/sum(Detroit[,3])))/.25

如果我提供了不必要的信息，我深表歉意。我真的很想利用邻接矩阵来计算这些行的每个“捆绑”的 CGI。

如果你碰巧知道我是怎么开始的，那就太好了。

对于任何新手错误，我深表歉意，我是 R 新手！

编辑：

我已经想出了如何解决这个问题，但是为了清楚问题并回答评论中提出的一个问题，让我说一个多边形的邻域是它自己和它相邻的每个多边形的联合。在我上面给出的示例中，对于多边形“A”，这将是多边形“A”、“C”和“E”的并集

4

2 回答 2

2

目前尚不清楚您想如何利用邻接矩阵。

一个想法是将您的问题表述为图表。igraph 适用于操作相邻的边和顶点。

这是我的想法：

# I read the adjency matrix
POLYID.adjency <- read.table(text ='A B C D E  
A  0 0 1 0 1  
B  0 0 1 0 0    
C  1 1 0 1 0     
D  0 0 1 0 1     
E  1 0 0 1 0',header = TRUE)
# I create the graph
require(igraph)
g <- graph.adjacency(adjmatrix=POLYID.adjency)
V(g)$label <- V(g)$name

作为选项，您可以绘制它：

 plot(g)

在此处输入图像描述

现在我使用属性矩阵，为每条边创建一个属性（因为每一行都是一条边）

# 创建一个虚拟属性矩阵

POLYID.attributes <- read.table(text =' TOT L10K 10_15K 15_20K 
A 500   24     30     77
B 400   25     30     87
C 300   26     30     97
D 200   27     30     57
E 100   28     30     47',header = TRUE)

# I set the attributes
for(x in colnames(POLYID.attributes)){
   g <- set.vertex.attribute(g, name = x,
                         value=  POLYID.attributes[,x])

  }

现在所有问题信息都在图表中。

str(g)
IGRAPH DN-- 5 10 -- 
+ attr: name (v/c), label (v/c), TOT (v/n), L10K (v/n),
        X10_15K (v/n), X15_20K (v/n)
+ edges (vertex names):
 [1] A->C A->E B->C C->A C->B C->D D->C D->E E->A E->D

现在我可以使用 igraph 选项获取每个节点的信息，例如：

e. 获取与B相邻的多边形的L10K属性

V(g)[get.adjlist(g,'out')$B]$L10K
[1] 26

在这里，我计算与 plyogon A 相邻的所有多边形的总和：

 sum(V(g)[get.adjlist(g,'out')$A]$TOT)
  400
于 2012-12-14T05:01:03.277 回答
0

这就是我最终要做的，虽然它看起来不像@agstudy 的那么优雅：

for(k in 1:nrow(adjacency00))
 {
  positions <- grep(1,adjacency00[k,])-1
  nghbrd <- UNSD00[c(positions,k),]

等，从而创建一个相邻多边形的框架，在其上进行后续计算

于 2012-12-17T18:59:24.563 回答