0

我有研究站点，我收集了数据和附近的气象站，其中包含有关温度和降水的信息。我想将我在研究地点的日常数据与最近气象站的天气信息配对。我认为，要做到这一点，我需要一个两步过程，首先选择离研究地点最近的气象站，然后使用天气数据创建一个新变量。

这是我的数据的快照：

# study sites
site <- rep(LETTERS[1:3], 5)
siteLat <- rep(c(41, 42, 44), 5)
siteLon <- rep(c(68, 62, 63), 5)
siteDate <- rep(1:5, 3)
dfSites <- data.frame(cbind(site, siteLat, siteLon, siteDate))

# weather stations
station <- rep(letters[1:3], 5)
stationLat <- rep(c(40, 43, 45), 5)
stationLon <- rep(c(67, 61, 64), 5)
stationDate <- rep(1:5, 3)
temp <- sample(10:20, 15, replace=TRUE)
dfStation <- data.frame(cbind(station, stationLat, stationLon, stationDate, temp))

我正在尝试使用这条线来确定最近的车站，但我只得到一排距离。

distVincentyEllipsoid(df2[c("recvLon", "recvLat")], weather[c("lon", "lat")])

计算完所有距离后，我对接下来的步骤有点不确定，但我想我需要一些东西来选择最近的车站和匹配日期。这是我想出的最好的：

dfSites %>% 
    mutate(closestStation = ???,
           temp1 = temp[station == closestStation & stationDate == siteDate])

最终结果是我的研究站点数据框，其中包含最近气象站的温度附加列。

4

1 回答 1

0

我认为distVincentyEllipsoid(p1, p2, ...)试图找到 的第一个点p1与 的第一个点之间的距离p2，第二个p1与 的第二个p2等的距离。您需要的是沿着 *"first in p1 对 all of 的扩展p2，第二个p1与 all of p2， ETC）。

调整您的代码以调用dfSitesand dfStation（而不是df2/ weather），以下内容应该适合您。（我将删除其中一个站点，dfStation[-1,...] 只是为了清楚地识别哪个维度代表站点与站点。

alldists <- sapply(seq_len(nrow(dfSites)), function(i) {
  distVincentyEllipsoid(dfSites[i,c("siteLon","siteLat")],
                        dfStation[-1,c("stationLon","stationLat")])
})
alldists
#           [,1]     [,2]     [,3]     [,4]     [,5]     [,6]     [,7]     [,8]
#  [1,] 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1
#  [2,] 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4
#  [3,] 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7
#  [4,] 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1
#  [5,] 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4
#  [6,] 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7
#  [7,] 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1
#  [8,] 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4
#  [9,] 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7
# [10,] 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1
# [11,] 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4
# [12,] 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7
# [13,] 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1
# [14,] 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4
#           [,9]    [,10]    [,11]    [,12]    [,13]    [,14]    [,15]
#  [1,] 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0
#  [2,] 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2
#  [3,] 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5
#  [4,] 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0
#  [5,] 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2
#  [6,] 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5
#  [7,] 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0
#  [8,] 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2
#  [9,] 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5
# [10,] 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0
# [11,] 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2
# [12,] 484015.5 119427.7 565573.7 484015.5 119427.7 565573.7 484015.5
# [13,] 228960.0 786180.9 123505.1 228960.0 786180.9 123505.1 228960.0
# [14,] 122086.2 481351.6 269760.4 122086.2 481351.6 269760.4 122086.2

（因为我们有14行，每一行都是你的一个站。你不应该做索引[-1,]，只知道哪一行/列。）由此，我们知道站点A和站之间的差异b是481351.6米（第一列，第二行）。

从这里，只需找到列最小值：

apply(alldists, 2, which.min)
#  [1] 3 1 2 3 1 2 3 1 2 3 1 2 3 1 2

建议离站点最近的站点A是b（which.min将返回第一个最小值，它不表示关系）。

现在，dfStation[apply(alldists, 2, which.min),]为您提供 15 行站数据，可以轻松cbind编辑或以其他方式与dfSites.

dplyr选项：

dfSites %>%
  mutate(
    station_i = purrr::map2_int(
      siteLat, siteLon,
      ~ which.min(geosphere::distVincentyEllipsoid(
          cbind(.x,.y), dfStation[-1,c("stationLon","stationLat")]))
      ),
    station = as.character(dfStation$station)[ station_i ]
  )
#    site siteLat siteLon siteDate station_i station
# 1     A      41      68        1         3       c
# 2     B      42      62        2         1       a
# 3     C      44      63        3         2       b
# 4     A      41      68        4         3       c
# 5     B      42      62        5         1       a
# 6     C      44      63        1         2       b
# 7     A      41      68        2         3       c
# 8     B      42      62        3         1       a
# 9     C      44      63        4         2       b
# 10    A      41      68        5         3       c
# 11    B      42      62        1         1       a
# 12    C      44      63        2         2       b
# 13    A      41      68        3         3       c
# 14    B      42      62        4         1       a
# 15    C      44      63        5         2       b

通过对它们进行外积可以看到轻微的（10-15%）速度提高。

outer(seq_len(nrow(dfSites)), seq_len(nrow(dfStation)),
      function(i,j) geosphere::distVincentyEllipsoid(dfSites[i,2:3], dfStation[j,2:3]))

这还会返回一个mxn矩阵（站行），然后您可以apply(...)通过该矩阵获取最近的索引。（我希望获得更大的性能提升，因为distVincentyEllipsoid只调用一次......）

于 2019-07-16T17:51:18.700 回答