R - 最大化多个场景的曲线下面积

R - maximize area under curve for multiple scenarios

考虑到我有两个向量,一个叫residues,另一个叫scores,它们有31个分数,每个残基一个,都是正数。为了说明,得到两个向量如下图:

residues <- 1:31
scores <- runif(n = 31, min = 0.35, max = 3.54)

我正在考虑一个随机序列来举例说明。 如果我绘制 residues x scores 我将得到以下图形:

我想做的是:我将考虑 15 个残基的特定组合(以下称为 15mer),跳过一个残基(即 1:15、2:16、3:17 一直到 17:31),我想计算所有这 17 种组合的曲线下面积 (AUC)。我的最终目标是 select 具有最高 AUC 的 15mer。

可以使用 zoo 包中的 rollmean 函数计算 AUC,如 this question 所示。然而,在这个例子中,我有 17 种可能的组合,我试图找到一个脚本来自动化这个过程。 提前致谢。

下面是我用过的函数

scores <- runif(n = 31, min = 0.35, max = 3.54)

fun <- function(dat, n) {
  require(zoo)
  N <- which(max(rollmean(dat, n)) == rollmean(dat, n))
  output <- matrix(0, length(N), n)
  for (i in 1:length(N)) {
   output[i, ] <- dat[N[i]:(N[i] + n - 1)]
  }
  output
}

fun(scores, 15)

让运行彻底翻过来

rollmean(dat, n)

您提到的 zoo 包为我们提供了滚动平均值,我们

max(rollmean(dat, n))

找到滚动平均值的最大值

max(rollmean(dat, n)) == rollmean(dat, n)

returns 一个TRUE/FALSE向量,其滚动均值等于最大值

N <- which(max(rollmean(dat, n)) == rollmean(dat, n))

returns 最大值的索引。根据您的数据,您可能有多个序列获得最大值,我们决定 return 所有这些都使用以下循环

for (i in 1:length(N)) {
  output[i, ] <- dat[N[i]:(N[i] + n -1)]
}

结果:

set.seed(12345)
scores <- runif(n = 31, min = 0.35, max = 3.54)

fun(scores, 15)
         [,1]     [,2]      [,3]     [,4]     [,5]    [,6]
[1,] 1.588179 1.633928 0.9208938 3.385791 1.797393 1.39234
         [,7]     [,8]     [,9]    [,10]    [,11]    [,12]
[1,] 3.429675 2.606867 2.406091 1.593553 2.578354 2.085545
       [,13]    [,14]    [,15]
[1,] 1.07243 1.895739 2.879693

fun(rpois(1000, 1), 10)
    [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,]    1    1    4    2    1    1    3    3    2     2
[2,]    1    4    2    1    1    3    3    2    2     1
[3,]    4    2    1    1    3    3    2    2    1     1
library(zoo)

set.seed(555)
residues <- 1:31
scores <- runif(n = 31, min = 0.35, max = 3.54)


which.max(sapply(1:17, function(x){sum(diff(residues[x:(x+14)])*rollmean(scores[x:(x+14)],2))}))
# result 7 i.e. 7:21

sapply(1:17, function(x){sum(diff(residues[x:(x+14)])*rollmean(scores[x:(x+14)],2))}) # gives you the AUCs
# result [1] 28.52530 29.10203 28.52847 27.65325 27.19925 28.77782 29.29373 28.13133 28.23705 27.68724 25.75294 25.27226 25.44963 25.81201 25.49907 23.48632
        #[17] 22.45763

或使用自定义函数

f_AUC <- function(x, y, lngth){
  sapply(1:(length(x)-lngth+1), function(z) sum(diff(x[z:(z+lngth-1)])*rollmean(y[z:(z+lngth-1)],2)))
}

f_AUC(x=residues, y=scores, lngth=15)