1

我想应用一个函数,将矩阵返回到大型 data.table 对象的每一行(原始文件大约 30 GB,我有 80 GB 内存),并取回一个 data.table 对象。我想有效地做到这一点。我目前的方法如下:

my.function <- function(x){
    alnRanges<-cigarToIRanges(x[6]);
    alnStarts<-start(alnRanges)+as.numeric(x[4])-1;
    alnEnds<-end(alnRanges)+as.numeric(x[4])-1;
    y<-x[-4];
    ys<-matrix(rep(y,length(alnRanges)),nrow=length(alnRanges),ncol=length(y),byrow=TRUE);
    ys<-cbind(ys,alnStarts,alnEnds);
    return(ys);     # ys is a matrix
}

    my.dt<-fread(my.file.name);
    my.list.of.matrices<-apply(my.dt,1,my.function);
    new.df<-do.call(rbind.data.frame,my.list.of.matrices);
    colnames(new.df)[1:14]<-colnames(my.dt)[-4];
    new.dt<-as.data.table(new.df);

注意1:我指定 my.function 只是为了表明它返回一个矩阵,因此我的应用行是一个矩阵列表。

注意2:我不确定我正在执行的操作有多慢,但似乎我可以减少行数。例如,将数据框转换为大对象的数据表会很慢吗?

可重现的例子:

请注意,Arun 和 Roland 让我更加努力地思考这个问题,所以我仍在努力解决它......可能是我不需要这些线......

我想获取一个 sam 文件,然后创建一个新的坐标文件,其中每个读取都根据其 CIGAR 字段进行拆分。

My sam file:
qname   rname   pos cigar
2218    chr1    24613476    42M2S
2067    chr1    87221030    44M
2129    chr1    79702717    44M
2165    chr1    43113438    44M
2086    chr1    52155089    4M921N40M

code:

library("data.table");
library("GenomicRanges");

sam2bed <- function(x){
    alnRanges<-cigarToIRanges(x[4]);
    alnStarts<-start(alnRanges)+as.numeric(x[3])-1;
    alnEnds<-end(alnRanges)+as.numeric(x[3])-1;
    #y<-as.data.frame(x[,pos:=NULL]);
    #ys<-y[rep(seq_len(nrow(y)),length(alnRanges)),];
    y<-x[-3];
    ys<-matrix(rep(y,length(alnRanges)),nrow=length(alnRanges),ncol=length(y),byrow=TRUE);
    ys<-cbind(ys,alnStarts,alnEnds);
    return(ys);
}


sam.chr.dt<-fread(sam.parent.chr.file);
setnames(sam.chr.dt,old=c("V1","V2","V3","V4"),new=c("qname","rname","pos","cigar"));
bed.chr.lom<-apply(sam.chr.dt,1,sam2bed);
> bed.chr.lom
[[1]]
                           alnStarts  alnEnds   
[1,] "2218" "chr1" "42M2S" "24613476" "24613517"

[[2]]
                         alnStarts  alnEnds   
[1,] "2067" "chr1" "44M" "87221030" "87221073"

[[3]]
                         alnStarts  alnEnds   
[1,] "2129" "chr1" "44M" "79702717" "79702760"

[[4]]
                         alnStarts  alnEnds   
[1,] "2165" "chr1" "44M" "43113438" "43113481"

[[5]]
                               alnStarts  alnEnds   
[1,] "2086" "chr1" "4M921N40M" "52155089" "52155092"
[2,] "2086" "chr1" "4M921N40M" "52156014" "52156053"

bed.chr.df<-do.call(rbind.data.frame,bed.chr.lom);

> bed.chr.df
    V1   V2        V3 alnStarts  alnEnds
1 2218 chr1     42M2S  24613476 24613517
2 2067 chr1       44M  87221030 87221073
3 2129 chr1       44M  79702717 79702760
4 2165 chr1       44M  43113438 43113481
5 2086 chr1 4M921N40M  52155089 52155092
6 2086 chr1 4M921N40M  52156014 52156053

bed.chr.dt<-as.data.table(bed.chr.df);

> bed.chr.dt
     V1   V2        V3 alnStarts  alnEnds
1: 2218 chr1     42M2S  24613476 24613517
2: 2067 chr1       44M  87221030 87221073
3: 2129 chr1       44M  79702717 79702760
4: 2165 chr1       44M  43113438 43113481
5: 2086 chr1 4M921N40M  52155089 52155092
6: 2086 chr1 4M921N40M  52156014 52156053
4

1 回答 1

3

假设ff是你的data.table,这个怎么样?

splits <- cigarToIRangesListByAlignment(ff$cigar, ff$pos, reduce.ranges = TRUE)
widths <- width(attr(splits, 'partitioning'))
cbind(data.table(qname=rep.int(ff$qname, widths), 
            rname=rep.int(ff$rname, widths)), as.data.frame(splits))

   qname rname space    start      end width
1:  2218  chr1     1 24613476 24613517    42
2:  2067  chr1     2 87221030 87221073    44
3:  2129  chr1     3 79702717 79702760    44
4:  2165  chr1     4 43113438 43113481    44
5:  2086  chr1     5 52155089 52155092     4
6:  2086  chr1     5 52156014 52156053    40
于 2013-09-21T16:59:12.463 回答