我正在研究温度的时空观测值，并将其存储在大小为100 * 100 * 504（100 * 100网格，用于504个不同小时，代表21天）的数组中。我正在从这些观察值中计算不同的指标（在不同的时间段（3到21天）内，这显然需要一些时间，并且我正在寻找提高计算效率的方法。我不太习惯R，所以我不确定我在做什么是否是最有效的方法。
我想做的一件事是（对于每个单元）找到温度超过特定阈值的最长连续时间。这是我目前正在做的事情：
utci_test = array(runif(100*100*504, min = 18, max = 42), c(100,100,504))
to_hs = function(utci, period=1:length(utci[1,1,]), hs_threshold){
utci_hs = utci*0
utci_hs[which(utci > hs_threshold)] = 1
utci_hs[is.na(utci)] = 0
return(utci_hs)
}
max_duration_hs = function(utci_hs, period=1:length(utci_hs[1,1,]) ){
apply(utci_hs, MARGIN=c(1,2), FUN=function(x){
r = rle(x)
max(r$lengths[as.logical(r$values)], fill = 0)
})
}
查看所需的时间后，我注意到第二步需要花费一些时间（请记住，我总共必须重复执行此操作约8000次）
system.time(to_hs(utci_test, hs_threshold=32.0))
# utilisateur système écoulé
# 0.051 0.004 0.055
system.time(to_hs(utci_test, hs_threshold=32.0))
# utilisateur système écoulé
# 0.053 0.000 0.052
utci_test_sh = to_hs(utci_test, hs_threshold=32.0)
system.time(max_duration_hs(utci_test_sh))
# utilisateur système écoulé
# 0.456 0.012 0.468
因此，我想知道是否有更有效的方法来执行此操作，因为我认为转换为Rle对象可能效率不高？
通过编写自己的有效版本的rle()函数，可以有点麻烦，因为您知道要运行1，并且进行的比较少了。这使您的速度提高了大约2倍，在我的机器（通用的Macbook）上，平均时间缩短了大约250毫秒左右。
如果您必须执行8,000次，则可以通过并行处理在多核计算机上运行的代码来节省最多的时间，这在R中很容易做到（例如，检查parallel软件包）。
下面的代码为加速。
# generate data
set.seed(123)
utci_test <- array(runif(100*100*504, min = 18, max = 42), c(100,100,504))
# original functions
to_hs = function(utci, period=1:length(utci[1,1,]), hs_threshold){
utci_hs = utci*0
utci_hs[which(utci > hs_threshold)] = 1
utci_hs[is.na(utci)] = 0
return(utci_hs)
}
max_duration_hs = function(utci_hs, period=1:length(utci_hs[1,1,]) ){
apply(utci_hs, MARGIN=c(1,2), FUN=function(x){
r = rle(x)
max(r$lengths[as.logical(r$values)], fill = 0)
})
}
# helper func for rle
rle_max <- function(v) {
max(diff(c(0L, which(v==0), length(v)+1))) - 1
}
max_dur_hs_2 <- function(utci_hs) {
apply(utci_hs, MARGIN=c(1,2), FUN= rle_max)
}
# Check equivalence
utci_hs <- to_hs(utci = utci_test, hs_threshold = 32)
all.equal(max_dur_hs_2(utci_hs),
max_duration_hs(utci_hs))
#> [1] TRUE
# Test speed
library(microbenchmark)
microbenchmark(max_dur_hs_2(utci_hs),
max_duration_hs(utci_hs))
#> Unit: milliseconds
#> expr min lq mean median uq max
#> max_dur_hs_2(utci_hs) 216.1481 236.7825 250.9277 247.9918 262.4369 296.0146
#> max_duration_hs(utci_hs) 454.5740 476.5710 501.5119 489.9536 509.8750 774.9963
#> neval cld
#> 100 a
#> 100 b
由reprex软件包（v0.3.0）创建于2020-05-07
本文收集自互联网，转载请注明来源。
如有侵权，请联系 [email protected] 删除。
我来说两句