向量化包含循環和if子句的搜索函數

Question

我給出了兩個非常大的數據集，並且我一直在嘗試構建一個函數，該函數將從一個集合中找出某些關於其他數據的if子句的某些座標組。 我的問題是，我寫的功能很慢雖然我一直在閱讀某些問題類似的問題的答案，但我還沒有設法使其工作。
所以，如果我給出：

>head(CTSS)  
    V1  V2  V3 
1 chr1 564563 564598 
2 chr1 564620 564649 
3 chr1 565369 565404 
4 chr1 565463 565541 
5 chr1 565653 565697 
6 chr1 565861 565922 

和

> head(href) 
    chr  region start  end strand nu gene_id transcript_id 
1 chr1 start_codon 67000042 67000044  + . NM_032291  NM_032291 
2 chr1   CDS 67000042 67000051  + 0 NM_032291  NM_032291 
3 chr1  exon 66999825 67000051  + . NM_032291  NM_032291 
4 chr1   CDS 67091530 67091593  + 2 NM_032291  NM_032291 
5 chr1  exon 67091530 67091593  + . NM_032291  NM_032291 
6 chr1   CDS 67098753 67098777  + 1 NM_032291  NM_032291 

對於從HREF數據集我想找到的前兩個值在起始列每個值CTSS數據集的第3列小於或等於，並將其保留在新的數據框中。
環路我寫道：

y <- CTSS[order(-CTSS$V3), ]  
find_CTSS <- function(x, y) { 
    n <- length(x$start) 
    foo <- data.frame(matrix(0, n, 6)) 
    for (i in 1:n) 
    { 
     a <- which(y$V3 <= x$start[i]) 
     foo[i, ] = c(x$start[i], x$stop[i], y$V2[a[1]], y$V3[a[1]] , y$V2[a[2]], y$V3[a[2]]) 
    } 

print(foo) 

} 

Answer 1

您提供很少的數據（but see here），所以這是一個有點難以基準您的解決方案。看看下面的解決方案是否滿足您的需求。

#make some fake data 
href <- data.frame(start = runif(10), stop = runif(10), other_col = sample(letters, 10)) 
CTSS <- data.frame(col1 = runif(100), col2 = runif(100)) 

# for each row in href (but extract only stop and start columns) 
result <- apply(X = href[, c("start", "stop")], MARGIN = 1, FUN = function(x, ctss) { 
      criterion <- x["start"] #make a criterion 
      #see which values are smaller or equal to this criterion (and sort them) 
      extracted <- sort(ctss[ctss$col2 <= criterion, "col2"]) 
      #extract last and one to last value 
      get.values <- extracted[c(length(extracted) - 1, length(extracted))] 
      #put values in data frame 
      out <- as.data.frame(matrix(get.values, ncol = 2)) 
      return(out) 
     }, ctss = CTSS) 

#pancake a list into a data.frame 
result <- do.call("rbind", result) 

Answer 2

我不知道我會花多少時間專注於此，所以我會向前邁進。當這類問題在APL雜誌上收到單行答案時，我是一名APL傢伙。後來，我成爲了一名C++/STL人，並以新的着裝代碼學習了所有相同的東西。有時R會讓我認爲APL與PHP配合。

在這個問題中，數據框是一個分心。這是一個簡單的矢量搜索，有些粘貼在一起。

對於性能關鍵向量搜索，您需要findInterval。搜索範圍需要訂購。 search-fors可以以任意順序，但對於大型列表，您需要訂購。

V <- sort (runif(10*1000*1000)) 
    U <- sort (runif(10*1000*1000)) 
    W <- findInterval (U, V) 

這運行在三個羊羔尾巴的搖晃。現在你有一對整數。第一列是1:length(U)，第二列的值是W內的整數索引。

sum(u==sort(u)[sort.int (sort.int (u, index.return=TRUE)$ix, index.return=TRUE)$ix]) 

好的，我的APL腦幹有貢獻。答案是長度（u），並演示了「粘合在一起」所需的逆向排序。

令人興奮的事實：只有R中sort函數的特殊情況返回索引向量。在APL中，這是你從成績上升/成績下降的唯一答案。但是，嘿，這不像他們第一次做對了。

您必須修改findInterval的結果才能在匹配位置的小於一邊選擇兩個元素，並且必須撤消兩種類型才能粘合在一起。我懷疑你的運行時將會被這兩種類型（非常長的列表）所佔據，或者組裝最終的數據框（用於較小的列表）。在我的系統上，對長度爲100 * 1000 * 1000的數字列表進行排序開始出現問題。

夾在中間的findInterval的運行時間將是一片薄薄的生菜，這讓我想起爲什麼我不打算閒逛。

Answer 3

我看到你想要的主要是加速。借用Roman Lustrik的代碼，我看不出有什麼優勢可以在應用程序中進行排序。這真的會減慢速度。事實上，你想盡可能多地從apply（循環）中獲得。所以下面的運行速度要快得多。

#all code using Roman Lustrik's made up data 

CTSSs <- CTSS[order(CTSS$col2),] #presort CTSS 
result <- lapply(X = href$start, FUN = function(x, ctss) { 
    extracted <- ctss$col2[ctss$col2 <= x] 
    get.values <- tail(extracted,2) 
    out <- matrix(get.values, ncol = 2) 
    return(out)}, ctss = CTSSs) 
#pancake a list into a data.frame 
result <- as.data.frame(do.call("rbind", result)) 

或者，我可以進一步遵循向量化的精神，真正地獲得儘可能小的功能。

CTSSs <- CTSS[order(CTSS$col2),] #presort CTSS 
extracted <- lapply(href$start, function(x, ctss) { 
    ctss$col2[ctss$col2 <= x]}, ctss = CTSSs) 
get.values <- lapply(extracted, tail, n = 2) 
result <- t(sapply(get.values, matrix, ncol = 2)) 

#convert to a data.frame 
result <- as.data.frame(result) 

這可能會更快，或者你的情況也許沒有，但是，你應該需要添加一箇中間步驟，這也可能採取真正的矢量內置函數的優勢，說，如果你想在做數學在將它們放入數據框之前的值，然後您可以輕鬆地做到這一點。另外，你會注意到現在我可以在矩陣階段使用一個sapply/transpose，它比lapply/rbind更快。這通常是加速向量化你的代碼的地方，而不是僅僅圍繞它做一個大循環。 （順便說一句，它可以更容易地在你的思維的每一步檢查錯誤...也許這不是一個人呢。）

修訂：

在進一步的思考，我意識到這是可以完全矢量。下面的代碼將生成你想要比任何前面的例子快得多的東西。訣竅是使用cut（）和aggregate（）命令。

href <- href[order(href$start),] #just sorted so that the 0 at the beginning makes sense and the labels then match 
margin <- cut(CTSS$col2, breaks = c(0,href$start), labels = href$start, right = TRUE) 
result <- aggregate(col2 ~ margin, data = CTSS, FUN = function(x) tail(x,2)) 

可以重新格式化，結果如你所願得到你想要什麼，但應該做它的肉。您可能需要將margin列更改爲numeric，以便它與href $ start匹配，並在上面的中間示例中使用與sapply類似的代碼將上面的項目對列表轉換爲兩個單獨的列。這是循環或應用語句中的if（）語句，之前放慢了你的速度，而cut（）消除了這一點。