2

ここから htsjdk.jar によって提供されるメソッドにアクセスしようとしています: https://samtools.github.io/htsjdk/

ここに文書化されています： https://samtools.github.io/htsjdk/javadoc/htsjdk/index.html

ジソンを使用。バイナリ BAM ファイルの開始位置と終了位置を取得するには、BAM ファイル インデックス (BAI ファイル) にアクセス/クエリする方法が必要です。テスト BAM および BAI ファイルは、 https ://github.com/samtools/htsjdk/tree/master/testdata/htsjdk/samtools/BAMFileIndexTest から取得できます 。

Jython レジストリに配置した後の jython 2.7.0:

python.security.respectJavaAccessibility = false
#I did in the jython comandline:
import sys
sys.path.append("/usr/local/soft/picard_1.138/htsjdk-1.138.jar")
from htsjdk.samtools import *
from java.io import File
#the BAM index file + BAM files 
bai_fh = File("./index_test.bam.bai")
mydict = SAMSequenceDictionary()
bai_foo = DiskBasedBAMFileIndex(bai_fh, mydict)

bai_foo.getNumberOfReferences() などのいくつかのメソッドにアクセスできますが、目的のメソッド
getBinsOverlapping(int referenceIndex, int startPos, int endPos) は BrowseableBAMIndex インターフェイスにあります。

しかし、Jython で Java クラスをサブクラス化することになると、私は迷ってしまいます。このタスクは、特定のゲノム位置に対応する BAM ファイル チャンクのリストを取得することです。テスト BAM/BAI ファイル、つまり chrM 10000-15000 (染色体、start_pos、end_pos) の場合、htsjdk の代わりに市販の samtools スタンドアロン プログラムを使用して、11 のマッピングされた読み取りを取得します。

samtools view index_test.bam  chrM:10000-15000

助けてくれて本当にありがとうございます

ダレク

編集： Groovy部分を再Groovyバージョン：2.4.4

groovy -cp libs/htsjdk-1.138.jar test_htsjdk.groovy

#!/usr/bin/env groovy

import htsjdk.samtools.*

File bam_fh =  new File("./A.bam")
File bai_fh =  new File("./A.bam.bai")

def mydict = new SAMSequenceDictionary()
def bai_foo = new DiskBasedBAMFileIndex(bai_fh, mydict)
println bai_foo.getNumberOfReferences()

上記のコードは groovy で動作します。私の問題は、このコードが機能しないことではなく、BAI ファイル形式を扱う Java クラスからメソッドにアクセスする適切な方法がわからないことです。htsjdk/src/java/htsjdk/samtools/*java ファイル (repo@github からの git clone) で AbstractBAMFileIndex を検索しましたが、何をする必要があるかはまだ明確ではありません。

4

2 に答える 2

3

彼らのgithubには 1 つの例があります。小さなセクションを添付しましたが、SamReader の使用方法の例が他にもあります。

    /**
     * Broken down
     */
    final SamReaderFactory factory =
            SamReaderFactory.makeDefault().enable(SamReaderFactory.Option.VALIDATE_CRC_CHECKSUMS).validationStringency(ValidationStringency.LENIENT);

    final SamInputResource resource = SamInputResource.of(new File("/my.bam")).index(new URL("http://broadinstitute.org/my.bam.bai"));

    final SamReader myReader = factory.open(resource);

    for (final SAMRecord samRecord : myReader) {
        System.err.print(samRecord);
    }

Groovy は、形式に関係なく多数のファイルを操作するのに非常に優れています。

new File('my/directory/with/many/files').eachFile{File readFile ->
    final SamInputResource resource = SamInputResource.of(readFile).index(new URL('IGotLazy.bai'))
    final SamReader myReader = ... etc

}
于 2015-09-09T22:25:03.950 に答える
2

参照インデックスに何を入力すればよいかわかりませんでしたが、この種のデータを扱ったことのない私にとっては、印刷物は十分に良さそうに見えました。

import sys
sys.path.append("/home/shackle/sam-tools/samtools-htsjdk-f650176/dist/htsjdk-1.138.jar")
from htsjdk.samtools import *
from java.io import File
#the BAM index file + BAM files 
indexFile = File("/home/shackle/index_test.bam.bai")
bamFile = File("/home/shackle/index_test.bam")
sfr = SAMFileReader(bamFile, indexFile)
sfr.enableIndexCaching(True)
bbi = sfr.getBrowseableIndex()
for i in range(0,256):
    print "i = " , i
    bl = bbi.getBinsOverlapping(i, 10000, 15000)
    count = 0
    for bin in bl:
        print "bin.getBinNumber() = " , bin.getBinNumber()
        count = count + 1
        print "count = ", count
于 2015-09-09T23:33:49.403 に答える