I'm kind of new to data.tables and I have a table containing DNA genomic coordinates like this:

   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

I wrote a custom function that I want to apply to each row of my circa 2 million rows table, it uses GenomicFeatures' mapToTranscripts to retrieve the transcript ID for the position and the CDS coordinate. I want to add them to my table in two new columns, like this:

   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

The function is the following:

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)
}

Then, I do:

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

And I get this error, indicating that the function is returning two lists of shorter length than the original table, instead of one new element per row:

    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 assumed that using the .I operator would apply the function on a row basis and return one value per row. I also made sure the function was not returning empty values using the if statement.

Then I tried this mock version of the function:

get_feature <- function(dt) {

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

}

And called it like this:

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

It makes a 1 row data table, instead of one the original length. So I concluded that my function, or maybe the way I'm calling it, is collapsing the elements of the resulting vector somehow. What am I doing wrong?

Update (with solution): I was able to figure out after deeper research around StackOverflow. As is explained in this answer and in ?data.table, .I is only intended for use in j (as in DT[i,j,by=]). Therefore, by=.I is equivalent to by=NULL and the proper syntax is by=1:nrow(dt) in order to group by row number and apply the function row-wise.

Unfortunately, for my particular case, this is utterly inefficient and I calculated an execution time of 20 seconds for 100 rows. For my 36 million row dataset that takes 3 months to complete.

In my case, I had to give up and use the mapToTranscripts function on the entire table like this, which takes a couple of seconds and was obviously the intended use.

    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)]

Then map back to the genome using mapFromTranscripts from GenomicFeatures package so I could use a data.tables join to retrieve information from the original table, which was the intended purpose of what I was trying to do.