gIntersectを使用して、複数のSpatialPolygonsのペアごとの交差を検索しますか？

Question

私は、空間ポリゴンのリストのすべてのペアごとの交点を含む空間オブジェクトを生成したいと思います。 gIntersectを手動で2つのレイヤーの交点に使用するのは簡単ですが、すべてのペアワイズ交点を同時に検索したいと思います。私は13のポリゴンを持っているので、重なりをチェックするために156の異なる組み合わせがあります。 lapplyかforループのどちらかが動作するかもしれないようですが、空間ポリゴンのすべての可能な組み合わせのマトリックスが必要になると思います。ここで

は、データのサブサンプルのための要点である：

# Overlap in spring and summer of 2016 
library(dplyr) 
library(gdata) 
library(sp) 
library(rgdal) 
library(rgeos) 
library(raster) 
library(spatstat) 

#Import Data (from gist) 
df16<-df16.sub 
#Extract only locations between April 1 to August 1, 2016 
sub16<-df16[df16$loctime>"2016-03-31"&df16$loctime<"2016-08-01",] 
#Subset by age 
colts.16<-sub16[sub16$age=="colt",] 
df<-colts.16 
df$id<-as.factor(df$id) 

#extract coordinates 
coordinates(df)<-c("location.long", "location.lat") 
#define coordinate ref system 
proj4string(df)<-CRS("+proj=longlat +datum=WGS84") 
#reproject 
df<-spTransform(df, CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5 
         +lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs ")) 
colnames(df@coords)<-c("x", "y") 
sp.df<-as.data.frame(df) 
crane.names<-as.list(unique(sp.df$id)) 
testdat<-sp.df 

temp.crane<-NULL 
temp.crane.day<-NULL 
dat.summary<-NULL 
tempdat<-NULL 
firstday<-NULL 
lastday<-NULL 
dat.summary<-NULL 
i<-NULL 
j<-NULL 
#loop to condense points into just roost sites 
for(i in 1:length(unique(testdat$id))){ 
    temp.crane<-testdat[testdat$id==crane.names[[i]],] 
    if(nrow(temp.crane)>0){ 
    current.crane<-as.character(crane.names[[i]]) 
    temp.days<-as.list(unique(temp.crane$day)) 
    days.vec<-unlist(temp.days) 
    firstday<-days.vec[1] 
    lastday<-last(days.vec) 
    for(j in seq.int(from=firstday, to=lastday)){ 
     tempdat<-temp.crane[which(temp.crane$day==j),]  
     firstrow<-tempdat[1,] 
     lastrow<-tempdat[nrow(tempdat),] 
     daylocs<-rbind(firstrow, lastrow) 

     if((j-1)<firstday) 
     {daylocs$night.dist<-NA} 
     else{ 
     prevdat<-temp.crane[which(temp.crane$day==j-1),]  
     firstrow<-prevdat[1,] 
     lastrow<-prevdat[nrow(prevdat),] 
     prevdaylocs<-rbind(firstrow, lastrow) 

     daylocs$night.dist<-sqrt((daylocs$x[1]-prevdaylocs$x[2])^2+(daylocs$y[1]-prevdaylocs$y[2])^2) 
     } 
      dat.summary<-rbind(dat.summary, daylocs) 
      } 
    } 
} 

roosts<-dat.summary 
temp.buffer<-NULL 
temp.crane<-NULL 
current.crane<-NULL 
temp.id<-NULL 
temp.sp.df<-NULL 
buff.spdf<-NULL 
crane.names<-as.list(unique(roosts$id)) 
test<-NULL 
#create list of SpatialPolygons 
for(i in 1:length(unique(roosts$id))){ 
    temp.crane<-roosts[roosts$id==crane.names[[i]],] 
    current.crane<-as.character(crane.names[[i]]) 
    coordinates(temp.crane)<-c("x","y") 
    proj4string(temp.crane)<-CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5 
          +lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs ") 
    temp.buffer<-gBuffer(temp.crane, width=3000) 
    test[i]<-list(temp.buffer) 
} 

私はその何かについての質問に適応しようとした。ここでは

https://gist.github.com/dwwolfson/b1dc7b9c084233a4a36401f7e7061897

は、私がこれまで試したものですコンパウンド（）関数を使用してちょっと関連しています。ここで、testは空間ポリゴンのリストです。

私は、次のエラーメッセージが出ます：次のエラーを生成しませんが、ArcMapの目視検査に、それだけで最後のラウンドはおそらく、重複のすべての機会を得ていない

Error in identical(spgeom1@proj4string, spgeom2@proj4string) :       
trying to get slot "proj4string" from an object of a basic class ("list") with no slots 

を。

overlap<-list() 
for(i in length(test)){ 
for(j in length(test)){ 
overlap<-gIntersection(test[[i]],test[[j]]) 
} 
} 

私の問題は、各反復のすべての空間ポリゴンオブジェクトを1つのリストに「rbind」またはマージする必要があることだと思います。 SpatialPolygonsの代わりにSpatialPolygonDataframesを最初に作成すると、組み合わせが簡単になるでしょうか？私は同じ番号であるiとjを避ける必要があります。そうでないと、同じポリゴンが重なり合うことになります。

次のコードを使用して、各繰り返しの重なりを空間オブジェクトのリストに格納しようとしましたが、結果を格納するオブジェクトがNULLとして開始されるため、最初の実行ではspRbindが異常終了します。

library(maptools) 
over.list<-NULL 

for(i in length(test)){ 
for(j in length(test)){ 
overlap<-gIntersection(test[[i]],test[[j]]) 
over.list<-spRbind(over.list, overlap) 
} 
} 

Answer 1

testの作成方法は非常に特殊です。これはSpatialオブジェクトではありません。
combnの使い方は正しくありません。このようなオブジェクトのリストは使用できません。ドキュメントは言う：

組み合わせ、または整数用のn X <ため
Xベクトル源combn - seq_len（n）を？。

ですので、そのまま使用してください。 あなたが与える例とデータはそのまま使用するのは難しいです。私は異なるものを修正しなければならなかった。次回は、データの問題に対応する簡単な例を提供する必要があります。参照先：How to make a great R reproducible example?
あなたの見ているものから理解できる限り、私はあなたに可能な方法を提案します。しかし、あなたはあなたのスクリプトをよりきれいにするために、R内のSpatialPolygonsのオブジェクトの構成が何であるかを知る必要があります。

ここでは、combnを使用しています。私はあなたがしたように、リスト内のすべてのポリゴンの交点を取得します。交差点がNULLでない場合は、私はまた、SpatialPolygonsとしてすべてのポリゴンを集めることができるようにするPolygonsオブジェクトとして取得します

combos <- combn(length(test),2) 
int <- vector("list", length = ncol(combos)) 
Polygons.list <- list() 
k.poly <- 0 # Number of non-null polygons 
for(k in 1:ncol(combos)){ 
    i<-combos[1,k] 
    j<-combos[2,k] 
    print(paste("intersecting", i, j)) 
    int.tmp <- gIntersection(test[[i]],test[[j]], byid=FALSE) 
    if (!is.null(int.tmp)) { 
    k.poly <- k.poly + 1 
    int[[k]] <- int.tmp 
    Polygons.list[[k.poly]] <- int[[k]]@polygons[[1]] 
    Polygons.list[[k.poly]]@ID <- paste(k.poly) 
    } 
} 

all.Intersections <- SpatialPolygons(Polygons.list) 
plot(all.Intersections, col = "red")