将自定义函数应用于 data.table 行 returns 不正确的值数量

Question

我是 data.tables 的新手，我有一个 table 包含这样的 DNA 基因组坐标：

       chrom   pause strand coverage
    1:     1 3025794      +        1
    2:     1 3102057      +        2
    3:     1 3102058      +        2
    4:     1 3102078      +        1
    5:     1 3108840      -        1
    6:     1 3133041      +        1

我写了一个自定义函数，我想将其应用于大约 200 万行的每一行 table，它使用 GenomicFeatures 的 mapToTranscripts 以字符串和新坐标的形式检索两个相关值。我想将它们添加到我的 table 的两个新列中，如下所示：

       chrom   pause strand coverage       transcriptID CDS
    1:     1 3025794      +        1 ENSMUST00000116652 196
    2:     1 3102057      +        2 ENSMUST00000116652  35
    3:     1 3102058      +        2 ENSMUST00000156816 888
    4:     1 3102078      +        1 ENSMUST00000156816 883
    5:     1 3108840      -        1 ENSMUST00000156816 882
    6:     1 3133041      +        1 ENSMUST00000156816 880

函数如下：

    get_feature <- function(dt){

      coordinate <- GRanges(dt$chrom, IRanges(dt$pause, width = 1), dt$strand) 
      hit <- mapToTranscripts(coordinate, cds_canonical, ignore.strand = FALSE) 
      tx_id <- tx_names[as.character(seqnames(hit))] 
      cds_coordinate <- sapply(ranges(hit), '[[', 1)

      if(length(tx_id) == 0 || length(cds_coordinate) == 0) {  
        out <- list('NaN', 0)
      } else {
        out <- list(tx_id, cds_coordinate)
      }

      return(out)
    }

那么，我做：

    counts[, c("transcriptID", "CDS"):=get_feature(.SD), by = .I]

我得到这个错误，表明该函数是 returning 两个长度比原始 table 更短的列表，而不是每行一个新元素：

Warning messages:
    1: In `[.data.table`(counts, , `:=`(c("transcriptID", "CDS"),  ... :
      Supplied 1112452 items to be assigned to 1886614 items of column 'transcriptID' (recycled leaving remainder of 774162 items).
    2: In `[.data.table`(counts, , `:=`(c("transcriptID", "CDS"),  ... :
      Supplied 1112452 items to be assigned to 1886614 items of column 'CDS' (recycled leaving remainder of 774162 items).

我假设使用 .I 运算符会逐行应用该函数，并且 return 每行一个值。我还使用 if 语句确保函数没有 returning 空值。

然后我尝试了这个函数的模拟版本：

    get_feature <- function(dt) {

      return('I should be returned once for each row')

    }

并这样称呼它：

    new.table <- counts[, get_feature(.SD), by = .I]

它生成 1 行数据 table，而不是原始长度。所以我得出结论，我的函数，或者我调用它的方式，正在以某种方式折叠结果向量的元素。我做错了什么？

更新（含解决方案）：正如@StatLearner指出的那样，在中解释说，如?data.table中所解释的那样，.I 仅适用于 j（如 DT[i,j,by=]）。因此，by=.I 等同于 by=NULL，正确的语法是 by=1:nrow(dt)，以便按行号分组并按行应用函数。

不幸的是，对于我的特殊情况，这是完全低效的，我计算出 100 行的执行时间为 20 秒。对于我需要 3 个月才能完成的 3600 万行数据集。

在我的例子中，我不得不放弃并在整个 table 上使用 mapToTranscripts 函数，这需要几秒钟，显然是预期的用途。

    get_features <- function(dt){
      coordinate <- GRanges(dt$chrom, IRanges(dt$pause, width = 1), dt$strand) # define coordinate
      hits <- mapToTranscripts(coordinate, cds_canonical, ignore.strand = FALSE) # map it to a transcript
      tx_hit <- as.character(seqnames(hits)) # get transcript number
      tx_id <- tx_names[tx_hit] # get transcript name from translation table

      return(data.table('transcriptID'= tx_id, 
                       'CDS_coordinate' =  start(hits))
    }

     density <- counts[, get_features(.SD)]

然后使用 GenomicFeatures 包中的 mapFromTranscripts 映射回基因组，这样我就可以使用 data.tables 连接从原始 table 中检索信息，这是我尝试做的事情的预期目的。

Answer 1

当我需要为 data.table 中的每一行应用一个函数时，我的做法是按行号分组：

counts[, get_feature(.SD), by = 1:nrow(counts)]

如中所述，.I 不适用于 by，因为它应该 return 由分组生成的行索引序列。 by = .I 不抛出错误的原因是 data.table 在 data.table 命名空间中创建对象 .I 等于 NULL，因此 by = .I 是等效的至 by = NULL.

请注意，使用 by=1:nrow(dt) 按行号分组并允许您的函数仅访问 data.table:

中的一行

require(data.table)
counts <- data.table(chrom = sample.int(10, size = 100, replace = TRUE),
                     pause = sample((3 * 10^6):(3.2 * 10^6), size = 100), 
                     strand = sample(c('-','+'), size = 100, replace = TRUE),
                     coverage = sample.int(3, size = 100, replace = TRUE))

get_feature <- function(dt){
    coordinate <- data.frame(dt$chrom, dt$pause, dt$strand)
    rowNum <- nrow(coordinate)
    return(list(text = 'Number of rows in dt', rowNum = rowNum))  
}

counts[, get_feature(.SD), by = 1:nrow(counts)]

将生成与 counts 中行数相同的 data.table，但 coordinate 将仅包含来自 counts

的一行

   nrow                 text rowNum
1:    1 Number of rows in dt      1
2:    2 Number of rows in dt      1
3:    3 Number of rows in dt      1
4:    4 Number of rows in dt      1
5:    5 Number of rows in dt      1

而 by = NULL 将向函数提供整个 data.table：

counts[, get_feature(.SD), by = NULL]

                   text rowNum
1: Number of rows in dt    100

这是 by 工作的预期方式。

将自定义函数应用于 data.table 行 returns 不正确的值数量

Applying custom function to data.table by row returns incorrect amount of values

r

data.table

genomicranges