부분 집합 SAM/BAM 파일 (R)

Question

많은 BAM 파일이 있습니다. 에서 scanBam으로 R에로드 할 수 있습니다.

그러나 읽기 전용 하위 집합 만 있으면됩니다. 나는에 관심이 있어요.

scanBam이의 읽기 모든 수천 데이터를 포함 13 개 요소 목록 1 개 요소 목록을 반환의 QName와 character 벡터를 가지고있다.

구조를 유지하는 qname으로이 개체의 하위 집합을 어떻게 만들 수 있습니까? 설명서 또는 온라인에서 아무것도 찾을 수 없습니다.

Answer 1

param = ScanBamParam (what = "qname")을 인수로 지정하여 qname을 포함하여 GenomicAlignments :: readGAlignments를 사용하여 데이터를 입력하는 것이 더 편리 할 것입니다. 그런 다음 % in %로 하위 집합을 구성 할 수 있습니다. 여기

library(GenomicAlignments) 
library(RNAseqData.HNRNPC.bam.chr14) 

fname <- RNAseqData.HNRNPC.bam.chr14_BAMFILES[1]  
want <- c("ERR127306.11930623", "ERR127306.24720935", 
    "ERR127306.23706011", "ERR127306.22418829", "ERR127306.13372247", 
    "ERR127306.20686334", "ERR127306.11412145", "ERR127306.4711647", 
    "ERR127306.7479582", "ERR127306.12737243") 
aln <- readGAlignments(fname, param=ScanBamParam(what="qname")) 
aln[mcols(aln)$qname %in% want] 

BAM 파일은 물론 큰, 그리고 QName의 그 큰 부분 인 ExperimentData 패키지 중 하나를 사용하여,보다 완벽한 예입니다; 파일을 반복하여 청크로 반복하는 것이 좋습니다. 이는 yieldReduce를 사용하여 (현재 Rsamtools에서) BamFile에 yieldSize를 합리적인 (예 : 1M) 읽기 횟수로 설정하고 MAP 함수로 데이터 청크를 입력하고 처리합니다 (예 : 원치 않는 읽기 필터링), 결과를 연결하는 (선택 사항) REDUCE 함수 및 반복이 완료된 시점을 나타내는 DONE (선택 사항) 함수가 있습니다. (yieldSize 샘플 데이터 그림 수 있도록 인위적으로 작고)와 같은 솔루션은 같습니다

bfl <- BamFile(fname, yieldSize=100000) ## larger, e.g., 1M-5M 
MAP <- function(bfl, want) { 
    ## message("tick") 
    aln <- readGAlignments(bfl, param=ScanBamParam(what="qname")) 
    if (length(aln) == 0) 
     NA       # semaphore -- DONE 
    else 
     aln[mcols(aln)$qname %in% want] 
} 
REDUCE <- c 
DONE <- function(x) identical(x, NA) 
result <- yieldReduce(bfl, MAP, REDUCE, DONE, want=want) 

하나는 scanBam를 사용하여 유사한 접근 방식을 채택 할 수있다,하지만 데이터 구조 (리스트의-목록)에 더 복잡하다 처리 대상 :

x <- scanBam(fname, param=ScanBamParam(what=c("qname", "pos"))) 
keep <- lapply(lapply(x, "[[", "qname"), "%in%", want) 
result <- Map(function(elts, keep) { 
    lapply(elts, "[", keep) 
}, x, keep) 

yieldReduce와 함께 사용할 수도 있습니다.

는 필터링과 새로운 빵 파일을 만드는 데 관심이 있다면

내가 subset(DataFrame(scanBam(bam_file)[[1]]),qname %in% qname_vector)를 사용하여 종료

filter_factory <- function(want) { 
    list(KeepQname = function(x) x$qname %in% want) 
} 
filter <- FilterRules(filter_factory(want)) 
dest <- filterBam(fname, tempfile(), filter=filter, 
        param=ScanBamParam(what="qname")) 
readGAlignments(dest) 

Answer 2

다음 읽습니다. 이것은 똑같은 구조 (목록의 목록)를 유지하지는 않지만 모든 정보가 유지되고 쉽게 액세스 할 수 있습니다.