Calcul de la distance {minimale} entre les polygones dans R

9

J'ai calculé la surface des distributions d'espèces (fusion de polygones à partir de fichiers de formes), mais comme cette zone peut être composée de polygones assez éloignés, je voudrais calculer une certaine mesure de dispersion. Ce que j'ai fait jusqu'à présent est de récupérer les centroïdes de chaque polygone, de calculer la distance entre ceux-ci et de les utiliser pour calculer le coefficient de variation, comme dans l'exemple factice ci-dessous;

require(sp)
require(ggplot2)
require(mapdata)
require(gridExtra)
require(scales)
require(rgeos)
require(spatstat)

# Create the coordinates for 3 squares
ls.coords <- list()
ls.coords <- list()
ls.coords[[1]] <- c(15.7, 42.3, # a list of coordinates
                    16.7, 42.3,
                    16.7, 41.6,
                    15.7, 41.6,
                    15.7, 42.3)

ls.coords[[2]] <- ls.coords[[1]]+0.5 # use simple offset

ls.coords[[3]] <- c(13.8, 45.4, # a list of coordinates
                    15.6, 45.4,
                    15.6, 43.7,
                    13.8, 43.7,
                    13.8, 45.4)

# Prepare lists to receive the sp objects and data frames
ls.polys <- list()
ls.sp.polys <- list()

for (ii in seq_along(ls.coords)) {
   crs.args <- "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
   my.rows <- length(ls.coords[[ii]])/2
   # create matrix of pairs
   my.coords <- matrix(ls.coords[[ii]],nrow = my.rows,ncol = 2,byrow = TRUE)
   # now build sp objects from scratch...
   poly = Polygon(my.coords)
   # layer by layer...
   polys = Polygons(list(poly),1)
   spolys = SpatialPolygons(list(polys))
   # projection is important
   proj4string(spolys) <- crs.args
   # Now save sp objects for later use
   ls.sp.polys[[ii]] <- spolys
   # Then create data frames for ggplot()
   poly.df <- fortify(spolys)
   poly.df$id <- ii
   ls.polys[[ii]] <- poly.df
}

# Convert the list of polygons to a list of owins
w <- lapply(ls.sp.polys, as.owin)
# Calculate the centroids and get the output to a matrix
centroid <- lapply(w, centroid.owin)
centroid <- lapply(centroid, rbind)
centroid <- lapply(centroid, function(x) rbind(unlist(x)))
centroid <- do.call('rbind', centroid)

# Create a new df and use fortify for ggplot
centroid_df <- fortify(as.data.frame(centroid))
# Add a group column
centroid_df$V3 <- rownames(centroid_df)

ggplot(data = italy, aes(x = long, y = lat, group = group)) +
  geom_polygon(fill = "grey50") +
  # Constrain the scale to 'zoom in'
  coord_cartesian(xlim = c(13, 19), ylim = c(41, 46)) +
  geom_polygon(data = ls.polys[[1]], aes(x = long, y = lat, group = group), fill = alpha("red", 0.3)) +
  geom_polygon(data = ls.polys[[2]], aes(x = long, y = lat, group = group), fill = alpha("green", 0.3)) +
  geom_polygon(data = ls.polys[[3]], aes(x = long, y = lat, group = group), fill = alpha("lightblue", 0.8)) + 
  coord_equal() +
  # Plot the centroids
  geom_point(data=centroid_points, aes(x = V1, y = V2, group = V3))

# Calculate the centroid distances using spDists {sp}
centroid_dists <- spDists(x=centroid, y=centroid, longlat=TRUE)

centroid_dists

       [,1]      [,2]     [,3]
[1,]   0.00000  69.16756 313.2383
[2,]  69.16756   0.00000 283.7120
[3,] 313.23834 283.71202   0.0000

# Calculate the coefficient of variation as a measure of polygon dispersion 
cv <- sd(centroid_dist)/mean(centroid_dist)
[1] 0.9835782

Tracé des trois polygones et de leurs centroïdes

entrez la description de l'image ici

Je ne sais pas si cette approche est très utile car dans de nombreux cas, certains des polygones (comme le bleu dans l'exemple ci-dessus) sont assez grands par rapport aux autres, augmentant ainsi encore plus la distance. Par exemple, le centroïde de l'Australie a presque la même distance à ses frontières occidentales qu'à Papau.

Ce que j'aimerais obtenir, c'est un apport sur des approches alternatives. Par exemple, comment ou avec quelle fonction puis-je calculer la distance entre les polygones?

J'ai testé pour convertir le cadre de données SpatialPolygon ci-dessus en PointPatterns (ppp) {spatstat}pour pouvoir exécuter nndist() {spatstat}pour calculer la distance entre tous les points. Mais étant donné que je traite de zones assez étendues (beaucoup de polygones et de grandes), la matrice devient énorme et je ne sais pas comment continuer à atteindre la distance minimale entre les polygones .

J'ai également regardé la fonction gDistance {rgeos}, mais je pense que cela ne fonctionne que sur des données projetées qui peuvent être un problème pour moi car mes zones peuvent en traverser plusieurs EPSG areas. Le même problème se poserait pour la fonction crossdist {spatstat}.

jO.
la source
1
Envisageriez-vous d'utiliser postgres/postgisen plus de R? J'ai utilisé un flux de travail dans lequel j'exécute la majorité de mon travail R, mais stocke les données dans une base de données à laquelle j'accède à l'aide sqldf. Cela vous permet d'utiliser toutes les postgisfonctions (dont la distance entre les polygones est simple)
djq
@djq: Merci d'avoir commenté. Ouais, j'essaierais certainement :) J'ai commencé à construire une base de données postgresmais j'ai arrêté quand je ne savais pas (ne regardais pas) comment connecter le workflow / géostats entre la base de données et R...
jO.

Réponses:

9

Vous pouvez faire cette analyse dans le package "spdep". Dans les fonctions voisines pertinentes, si vous utilisez "longlat = TRUE", la fonction calcule la distance du grand cercle et renvoie les kilomètres comme unité de distance. Dans l'exemple ci-dessous, vous pouvez contraindre l'objet de liste de distances résultant ("dist.list") à une matrice ou à data.frame, cependant, il est assez efficace de calculer des statistiques récapitulatives à l'aide de lapply.

require(sp)
require(spdep)

# Create SpatialPolygonsDataFrame for 3 squares
poly1 <- Polygons(list(Polygon(matrix(c(15.7,42.3,16.7,42.3,16.7,41.6,15.7,41.6,15.7,42.3), 
                   nrow=5, ncol=2, byrow=TRUE))),"1")     
poly2 <- Polygons(list(Polygon(matrix(c(15.7,42.3,16.7,42.3,16.7,41.6,15.7,41.6,15.7,42.3)+0.5, 
                   nrow=5, ncol=2, byrow=TRUE))),"2")     
poly3 <- Polygons(list(Polygon(matrix(c(13.8, 45.4, 15.6, 45.4,15.6, 43.7,13.8, 43.7,13.8, 45.4), 
                   nrow=5, ncol=2, byrow=TRUE))),"3")                      
spolys = SpatialPolygons(list(poly1,poly2,poly3),1:3)
 spolys <- SpatialPolygonsDataFrame(spolys, data.frame(ID=sapply(slot(spolys, "polygons"), 
                                    function(x) slot(x, "ID"))) )   
   proj4string(spolys) <- "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"

# Centroid coordinates (not used but provided for example) 
coords <- coordinates(spolys)

# Create K Nearest Neighbor list
skNN.nb <- knn2nb(knearneigh(coordinates(spolys), longlat=TRUE), 
                  row.names=spolys@data$ID)

# Calculate maximum distance for all linkages 
maxDist <- max(unlist(nbdists(skNN.nb, coordinates(spolys), longlat=TRUE)))

# Create spdep distance object
sDist <- dnearneigh(coordinates(spolys), 0, maxDist^2, row.names=spolys@data$ID)
  summary(sDist, coordinates(spolys), longlat=TRUE)

# Plot neighbor linkages                  
plot(spolys, border="grey") 
  plot(sDist, coordinates(spolys), add=TRUE)  

# Create neighbor distance list 
( dist.list <- nbdists(sDist, coordinates(spolys), longlat=TRUE) )

# Minimum distance 
( dist.min <- lapply(dist.list, FUN=min) )

# Distance coefficient of variation    
( dist.cv <- lapply(dist.list, FUN=function(x) { sd(x) / mean(x) } ) )
Jeffrey Evans
la source
Merci pour vos commentaires et la perspicacité dans le spdebpaquet. Juste pour clarifier, cette approche donne le même résultat que dans mon exemple, non?
jO.
Juste au cas où vous ne
verriez
Bien que la réponse fournisse un code utile pour calculer les distances entre les centroïdes, elle ne traite pas du point central du PO qui est de savoir comment trouver la distance entre les deux points les plus proches des frontières du polygone.
csfowler du
Un énorme flic et une mauvaise forme pour SE, mais je ne peux pas faire tout le travail en ce moment. Ma propre recherche de réponse à cette question semble indiquer que la fonction gDistance des bibliothèques de bibliothèque fera ce que l'OP voulait: trouver la distance la plus courte entre les bords. Si, dans ma hâte de respecter un délai serré, j'ai mal interprété OP ou Jeffrey Evans mes sincères excuses.
csfowler