Rleベクターから効率的にGRanges/IRangesを構築する

Question

私は、ゲノム上のすべての位置でいくつかの値を順番に表すランレングスコード化ベクトルを持っています。おもちゃの例として、私は長さ10のただ1つの染色体を持っていた、そして私はグランジオブジェクトにこれを強制したい

library(GenomicRanges) 

set.seed(1) 
toyData = Rle(sample(1:3,10,replace=TRUE)) 

ように見えるベクトルを持っているでしょうとします。私が思いつくことができる最高は

gr = GRanges('toyChr',IRanges(cumsum(c(0,runLength(toyData)[-nrun(toyData)])), 
           width=runLength(toyData)), 
      toyData = runValue(toyData)) 

ですが、動作は遅いです。同じオブジェクトをより速く構築する方法はありますか？

Answer 1

@TheUnfunCatが指摘したように、OPのソリューションはかなり安定しています。以下のソリューションは、元のソリューションよりも適度に高速です。私はbase Rのほぼすべての組み合わせを試して、効率RleをS4Vectorsパッケージから叩くことができなかったので、私はRcppに頼った。ここでは主な機能は次のとおりです。

GenomeRcpp <- function(v) { 
    x <- WhichDiffZero(v) 
    m <- v[c(1L,x+1L)] 
    s <- c(0L,x) 
    e <- c(x,length(v))-1L 
    GRanges('toyChr',IRanges(start = s, end = e), toyData = m) 
} 

WhichDiffZeroはかなりbase Rでwhich(diff(v) != 0)とまったく同じことを行いRcpp機能です。多くのクレジットは@G.Grothendieckになります。以下は

#include <Rcpp.h> 
using namespace Rcpp; 

// [[Rcpp::export]] 
IntegerVector WhichDiffZero(IntegerVector x) { 
    int nx = x.size()-1; 
    std::vector<int> y; 
    y.reserve(nx); 
    for(int i = 0; i < nx; i++) { 
     if (x[i] != x[i+1]) y.push_back(i+1); 
    } 
    return wrap(y); 
} 

いくつかのベンチマークです：

set.seed(437) 
testData <- do.call(c,lapply(1:10^5, function(x) rep(sample(1:50, 1), sample(1:30, 1)))) 

microbenchmark(GenomeRcpp(testData), GenomeOrig(testData)) 
Unit: milliseconds 
       expr  min  lq  mean median  uq  max neval cld 
GenomeRcpp(testData) 20.30118 22.45121 26.59644 24.62041 27.28459 198.9773 100 a 
GenomeOrig(testData) 25.11047 27.12811 31.73180 28.96914 32.16538 225.1727 100 a 

identical(GenomeRcpp(testData), GenomeOrig(testData)) 
[1] TRUE 

私は過去数日間、このオフにし、に取り組んできたと私は間違いなく満足していません。私は誰かが私がやったことを（別のアプローチであるので）取って、もっと良いものを作ることを望んでいます。