有很多種格式來記錄mapping的結果，大多數都收錄在UCSC的幫助文檔中。在mapping之后，如何對mapping的結果進行評估呢？ 最簡單的，就是通過samtools來評估mapping質量了。

samtools idxstats aln.sorted.bam

注意，這一步之前需要經過sort和index。結果會顯示：

chr1 195471971 6112404 0chr10 130694993 3933316 0chr11 122082543 6550325 0chr12 120129022 3876527 0chr13 120421639 5511799 0chr14 124902244 3949332 0chr15 104043685 3872649 0chr16 98207768 6038669 0chr17 94987271 13544866 0chr18 90702639 4739331 0chr19 61431566 2706779 0chr2 182113224 8517357 0chr3 160039680 5647950 0chr4 156508116 4880584 0chr5 151834684 6134814 0chr6 149736546 7955095 0chr7 145441459 5463859 0chr8 129401213 5216734 0chr9 124595110 7122219 0chrM 16299 1091260 0chrX 171031299 3248378 0chrY 91744698 259078 0* 0 0 0

其中第一列是染色體名稱，第二列是序列長度，第三列是mapped reads數，第四列是unmapped reads數。 如果是RNAseq，我們可以使用broad institute的RNA-SeQC來得到更加完整的報告。下載到文件之后，也許需要安裝BWA來獲取更精準的結果，但是如果不安裝的話，也可以進行分 析。一般來說，這一步不需要特別精準的結果，所以我很少使用BWA選項。下載的文件如果是.zip結尾的，直接把它改寫成.jar就可以運行了。 在它的主頁上下載所需要的Example RNA-seq Data。下載結束之后，該解壓的解壓縮。接下來運行：

samtools index example/ThousandReads.bamsamtools faidx example/Homo_sapiens_assembly19.fastajava -Xmx2048m -jar RNA-SeQC_v1.1.7.jar -n 1000 -s "TestId|example/ThousandReads.bam|TestDesc" -t example/gencode.v7.annotation_goodContig.gtf -r example/Homo_sapiens_assembly19.fasta -o ./testReport/ -start gc -gc example/gencode.v7.gc.txt

以上的參數只有一個與其說明文檔不一樣的地方就是使用了-Xmx2048m來指定java虛擬機的內存大小為2G。如果遇到java.lang.OutOfMemoryError，還可以指定得再大些。

當然如果是自己的文件的話，還需要多兩步：

1.BAM，reference及GTF文件的基因組名稱必須一致。

2.需要使用picard工具包中的CreateSequenceDictionary來構建一個dictionary文件。

轉載于:https://www.cnblogs.com/pennyy/p/4650690.html

總結

以上是生活随笔為你收集整理的samtools idxstats的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。