问题描述
我在R中有一组两个数据帧
First:
site_no <- c("02110500","02110550","02110701","02110704","02110760","02110777","021108044","02110815")
lat_coor <- c(33.91267,33.85083,33.86100,33.83295,33.74073,33.85156,33.65017,33.44461)
long_coor <- c(-78.71502,-78.89722,-79.04115,-79.04365,-78.86669,-78.65585,-79.12310,-79.17393)
AllStations <- data.frame(site_no,lat_coor,long_coor)
Second:
station <- c("USGS-02146110","USGS-02146110","USGS-021473426","USGS-021473426")
latitude <- c(34.88928,34.85651,34.71679,34.24320,34.80012,34.80012)
longitude <- c(-81.06869,-82.22622,-82.17372,-81.31954,-82.36512,-82.36512)
ContaminantStations <- data.frame(station,latitude,longitude)
我的数据集要长得多,但是出于这个问题的目的,我认为这应该足够了。
我想从第一个数据框(AllStations)中找到第二个数据框(ContaminantStations)中点半径内的所有站,并将它们附加到新的数据框中(仅来自AllStations),我需要提取包含所有信息的站点。我尝试了一些合乎逻辑的方法,但是它们都不起作用或没有道理。我也尝试使用RANN:nn2,但这只为我计数。
任何帮助将不胜感激
解决方法
我认为您只需要遍历AllStations
内的每个对象,并返回半径内ContaminantStations
的最接近值。
func <- function(stations,constations,radius = 250000) {
if (!NROW(stations) || !NROW(constations)) return()
if (length(radius) == 1 && NROW(constations) > 1) {
radius <- rep(radius,NROW(constations))
} else if (length(radius) != NROW(constations)) {
stop("'radius' must be length 1 or the same as the number of rows in 'constations'")
}
out <- integer(NROW(stations))
for (i in seq_len(NROW(stations))) {
dists <- geosphere::distHaversine(stations[i,],constations)
out[i] <- if (any(dists <= radius)) which.min(dists) else 0L
}
return(out)
}
这将返回一个整数向量,指示最接近的污染工位。如果半径范围内都不存在,则返回0
。安全地将其用作原始帧上的行索引。
每个参数只能包含两列,第一列是经度。 (我不对函数中的列名做任何假设。)radius
以米为单位,与geosphere
包假设一致。
ind <- func(AllStations[,c("long_coor","lat_coor")],ContaminantStations[,c("longitude","latitude")],radius = 230000)
ind
# [1] 0 6 6 6 0 0 6 6
这些是ContaminantStations
行上的索引,其中非零表示污染站最接近AllStations
的特定行。
我们可以确定与此最接近的污染物站(有很多方法可以做到这一点,包括tidyverse和其他技术……这只是一个开始)。
AllStations$ClosestContaminantStation <- NA_character_
AllStations$ClosestContaminantStation[ind > 0] <- ContaminantStations$station[ind]
AllStations
# site_no lat_coor long_coor ClosestContaminantStation
# 1 02110500 33.91267 -78.71502 <NA>
# 2 02110550 33.85083 -78.89722 USGS-021473426
# 3 02110701 33.86100 -79.04115 USGS-021473426
# 4 02110704 33.83295 -79.04365 USGS-021473426
# 5 02110760 33.74073 -78.86669 <NA>
# 6 02110777 33.85156 -78.65585 <NA>
# 7 021108044 33.65017 -79.12310 USGS-021473426
# 8 02110815 33.44461 -79.17393 USGS-021473426
查看数据的依据:
这种方法的另一种方法是返回最近的污染物测站的距离和索引,而与半径无关,以便稍后进行过滤。
func2 <- function(stations,NROW(constations))
} else if (length(radius) != NROW(constations)) {
stop("'radius' must be length 1 or the same as the number of rows in 'constations'")
}
out <- data.frame(ind = integer(NROW(stations)),dist = numeric(NROW(stations)))
for (i in seq_len(NROW(stations))) {
dists <- geosphere::distHaversine(stations[i,constations)
out$ind[i] <- which.min(dists)
out$dist[i] <- min(dists)
}
return(out)
}
示范,包括将污染站置于同一帧中。
AllStations2 <- cbind(
AllStations,func2(AllStations[,"latitude")])
)
AllStations2
# site_no lat_coor long_coor ind dist
# 1 02110500 33.91267 -78.71502 1 241971.5
# 2 02110550 33.85083 -78.89722 6 227650.6
# 3 02110701 33.86100 -79.04115 6 214397.8
# 4 02110704 33.83295 -79.04365 6 214847.7
# 5 02110760 33.74073 -78.86669 6 233190.8
# 6 02110777 33.85156 -78.65585 6 249519.7
# 7 021108044 33.65017 -79.12310 6 213299.3
# 8 02110815 33.44461 -79.17393 6 217378.9
AllStations3 <- cbind(
AllStations2,ContaminantStations[AllStations2$ind,]
)
AllStations3
# site_no lat_coor long_coor ind dist station latitude longitude
# 1 02110500 33.91267 -78.71502 1 241971.5 USGS-02146110 34.88928 -81.06869
# 6 02110550 33.85083 -78.89722 6 227650.6 USGS-021473426 34.24320 -81.31954
# 6.1 02110701 33.86100 -79.04115 6 214397.8 USGS-021473426 34.24320 -81.31954
# 6.2 02110704 33.83295 -79.04365 6 214847.7 USGS-021473426 34.24320 -81.31954
# 6.3 02110760 33.74073 -78.86669 6 233190.8 USGS-021473426 34.24320 -81.31954
# 6.4 02110777 33.85156 -78.65585 6 249519.7 USGS-021473426 34.24320 -81.31954
# 6.5 021108044 33.65017 -79.12310 6 213299.3 USGS-021473426 34.24320 -81.31954
# 6.6 02110815 33.44461 -79.17393 6 217378.9 USGS-021473426 34.24320 -81.31954
您可以从此处随意选择半径:
subset(AllStations3,dist < 230000)
# site_no lat_coor long_coor ind dist station latitude longitude
# 6 02110550 33.85083 -78.89722 6 227650.6 USGS-021473426 34.2432 -81.31954
# 6.1 02110701 33.86100 -79.04115 6 214397.8 USGS-021473426 34.2432 -81.31954
# 6.2 02110704 33.83295 -79.04365 6 214847.7 USGS-021473426 34.2432 -81.31954
# 6.5 021108044 33.65017 -79.12310 6 213299.3 USGS-021473426 34.2432 -81.31954
# 6.6 02110815 33.44461 -79.17393 6 217378.9 USGS-021473426 34.2432 -81.31954