微信公众号搜"智元新知"关注
微信扫一扫可直接关注哦!

对于地理点纬度/经度的全球数据集,如何找到最近的邻居,以解释我们星球的球形特性

如何解决对于地理点纬度/经度的全球数据集,如何找到最近的邻居,以解释我们星球的球形特性

所以我们在不考虑投影问题的情况下完成了这项工作。问题是在哪里(以及如何)最好地添加重新投影,以便函数返回以公里为单位的值而不是当前度数:

library(raster)
library(purrr)
library(sf)

#example presence data from model
r1 <- raster(nrow=360,ncol=720)
crs(r1) <- "+proj=longlat +datum=wgs84 +no_defs"
values(r1) <- rbinom(ncell(r1),2,0.01)

r1_points <- rasterToPoints(r1)
r1_df <- data.frame(r1_points)
r1_presence <- r1_df %>% dplyr::filter(layer==1)

#example survey data
survey_points <- cbind(rnorm(50) * 5 + 10,rnorm(50) + 50)
pt2 <- st_multipoint(cbind(survey_points[,1],survey_points[,2]))

#distance between each modelled presence (pt1) and survey point (pt2)
get_distances <- function(i,pt2,df) {
  pt1 <- st_multipoint(cbind(df[i,df[i,2]),dim = "XY")
  a <- st_nearest_points(pt1,pt2)
  return(st_length(a))
}

#loop for all modelled presences
output <- map_dbl(1:nrow(r1_presence),get_distances,r1_presence)

理想情况下,一个完美的答案是扩展 get_distances 函数添加一个新选项,该选项执行适当的重新投影并以公里为单位返回值。

这里可能有几种不同的方法,我很好奇人们会想出什么。

解决方法

你问题的前提是错误的。通常不应将数据投影到计算距离;投影会失真,因此距离不会很精确。如果您的数据具有较大的空间范围,则尤其如此。

计算最近距离的方法可能会有所不同,例如取决于您拥有的点数。但在这个例子中,你可以只使用蛮力并计算所有距离

library(raster)

r1 <- raster(nrow=360,ncol=720)
crs(r1) <- "+proj=longlat +datum=WGS84 +no_defs"
set.seed(1)
values(r1) <- rbinom(ncell(r1),2,0.01)
r1_presence <- rasterToPoints(r1,fun=function(x)x==1)
survey_points <- cbind(rnorm(50) * 5 + 10,rnorm(50) + 50)

d <- pointDistance(r1_presence[,1:2],survey_points,lonlat=TRUE)
mind <- apply(d,1,min)

head(mind)
#[1] 4268674 4261209 4258182 4254560 4220200 4218188

使用 terra v 1.2-0(目前是开发版本,您可以install from github)可以做到

library(terra)
#terra version 1.2.0
from <- vect(r1_presence[,crs="+proj=longlat +datum=WGS84")
to <- vect(survey_points,crs="+proj=longlat +datum=WGS84")
x <- nearest(from,to)
x
# class       : SpatVector 
# geometry    : points 
# dimensions  : 5158,7  (geometries,attributes)
# extent      : -2.271833,22.16701,48.51506,51.98087  (xmin,xmax,ymin,ymax)
# coord. ref. : +proj=longlat +datum=WGS84 +no_defs 
# names       : from_id from_x from_y to_id  to_x  to_y  distance
# type        :   <num>  <num>  <num> <num> <num> <num>     <num>
# values      :       1 -171.2  89.75    25 8.752 51.98 4.269e+06
#                     2 -128.2  89.75    25 8.752 51.98 4.261e+06
#                     3 -119.8  89.75    25 8.752 51.98 4.258e+06

版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。