R while循环与矢量条件

R while loop with vector condition

我想向量化一个使用 while 循环的函数。

原函数为

getParamsLeadtime <- function(leadtimeMean, in_tolerance, tolerance){
  searchShape=0
  quantil=0

  # iterates the parameters until the percentage of values is within the interval of tolerance
  while (quantil < in_tolerance){
    searchShape = searchShape+1
    quantil <- pgamma(leadtimeMean+tolerance,shape=searchShape,rate=searchShape/leadtimeMean) -
                  pgamma(leadtimeMean-tolerance,shape=searchShape,rate=searchShape/leadtimeMean) 
  }

  leadtimeShape <- searchShape
  leadtimeRate <- searchShape/leadtimeMean 

  return(c(leadtimeShape, leadtimeRate))
}

我想对该函数进行矢量化调用,以将其应用于数据框。目前我正在遍历它:

leadtimes <- data.frame()

for (a in seq(92:103)) {
  leadtimes <- rbind(leadtimes, getParamsLeadtime(a, .85,2))
}

当我尝试对函数进行向量化时,while 似乎不接受向量作为条件。出现以下警告:

Warning message:
In while (input["U"] < rep(tolerance, dim(input)[1])) { :
  the condition has length > 1 and only the first element will be used

这让我假设 while 不喜欢矢量。你能告诉我如何向量化函数吗?

在旁注中,我想知道为什么生成的提前期 -data.frame 的列名称看起来是值:

> leadtimes
   X1     X1.1
1   1 1.000000
2   1 0.500000
3   4 1.333333
4   8 2.000000
5  13 2.600000
6  19 3.166667
7  25 3.571429
8  33 4.125000
9  42 4.666667
10 52 5.200000
11 63 5.727273
12 74 6.166667

这是一个非常高效的选项。

对于给定的平均提前期,对于 +tol-tol 的情况,我们在足够大的 shp 序列上对 pgamma 的计算进行矢量化.我们计算(矢量化)差异,并与 in_tol 进行比较。大于 in_tol 的向量的第一个元素的索引(负 1,因为我们的序列从 0 开始)是导致 pgamma 大于 [=] 的 shp 的最低值16=].

f <- function(lead, in_tol, tol) {
  shp <- which(!(pgamma(lead + tol, 0:10000, (0:10000)/lead) - 
                 pgamma(lead - tol, 0:10000, (0:10000)/lead)) 
               < in_tol)[1] - 1
  rate <- shp/lead
  c(shp, rate)
}

然后我们可以 sapply 在平均交货时间范围内。

t(sapply(1:12, f, 0.85, 2))

##       [,1]     [,2]
##  [1,]    1 1.000000
##  [2,]    1 0.500000
##  [3,]    4 1.333333
##  [4,]    8 2.000000
##  [5,]   13 2.600000
##  [6,]   19 3.166667
##  [7,]   25 3.571429
##  [8,]   33 4.125000
##  [9,]   42 4.666667
## [10,]   52 5.200000
## [11,]   63 5.727273
## [12,]   74 6.166667

system.time(leadtimes <- sapply(1:103, f, 0.85, 2))

##   user  system elapsed 
##   1.28    0.00    1.30 

您只需要确保为形状参数选择一个合理的上限(这里我选择了 10000,这非常慷慨)。请注意,如果您 没有 选择足够高的上限,一些 return 值将是 NA.