How "bases mapped (cigar)" value in samtools stats is calculated?
Hi:
I am trying to use samtools stats to estimate "bases mapped (cigar)" value from some target regions. But I can't figure out how the number is calcualted out. Below is a small test:
This is a sam file subset "test_subset.sam", contains only 15 reads that all covers posiotion chr22:31656041. All 15 reads are 34M in cigar.:
@HD VN:1.3
@SQ SN:chr1 LN:249250621
@SQ SN:chr2 LN:243199373
@SQ SN:chr3 LN:198022430
@SQ SN:chr4 LN:191154276
@SQ SN:chr5 LN:180915260
@SQ SN:chr6 LN:171115067
@SQ SN:chr7 LN:159138663
@SQ SN:chr8 LN:146364022
@SQ SN:chr9 LN:141213431
@SQ SN:chr10 LN:135534747
@SQ SN:chr11 LN:135006516
@SQ SN:chr12 LN:133851895
@SQ SN:chr13 LN:115169878
@SQ SN:chr14 LN:107349540
@SQ SN:chr15 LN:102531392
@SQ SN:chr16 LN:90354753
@SQ SN:chr17 LN:81195210
@SQ SN:chr18 LN:78077248
@SQ SN:chr19 LN:59128983
@SQ SN:chr20 LN:63025520
@SQ SN:chr21 LN:48129895
@SQ SN:chr22 LN:51304566
@SQ SN:chrX LN:155270560
@SQ SN:chrM LN:16571
USI-EAS28:3:8:1631:168#0 16 chr22 31656010 37 34M * 0 0 CGTGTCCGTGGAGAGTGCCTGCTCCAACTACGCC \ZZNON\QTUYEGa``R^^``\^`a_`a`X^\`_ XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
USI-EAS28:1:82:1303:1490#0 16 chr22 31656011 37 34M * 0 0 GTGTCGGTGGAGAGTGCCTGCTCCAACTACGCCA BBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB XT:A:U NM:i:1 X0:i:1 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:5C28
USI-EAS28:2:36:650:781#0 16 chr22 31656012 37 34M * 0 0 TGTCCGTGGAGATTGCCTGCTCCAACTACGCCAC BBBS]NUOWRKLPSDYUXXXZVVVM\ZN___ZT\ XT:A:U NM:i:1 X0:i:1 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:12G21
USI-EAS28:3:71:1674:707#0 16 chr22 31656012 37 34M * 0 0 TGTCCGTGGAGAGTGCCTGCTCCAACTACGCCAC _]^]WY[a```a``\[UaaX`_aa`[`^aa`_a` XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
HWI-EAS276:4:52:1678:656#0 16 chr22 31656012 37 34M * 0 0 TGTCCGTGGAGAGTGCCTGCTCCAACTACGCCAC Z_```_U_a\a\aZa`\P__O_`]S]JV`a``Y` XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
HWI-EAS276:4:55:1651:892#0 16 chr22 31656012 37 34M * 0 0 TGTCCGTGGAGAGTGCCTGCTCCAACTACGCCAC [a\`__V\aTa_aaa`YL_aa^``WaR\aaaaZ_ XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
HWI-EAS276:4:57:1535:1224#0 16 chr22 31656014 37 34M * 0 0 TCCGTGGAGAGTGCCTGCTCCAACTACGCCACCG `aa_X_aaa`_Z_[\O]][[aa]a`\aaa_\`aa XT:A:U NM:i:1 X0:i:1 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:33A0
HWI-EAS276:4:6:1768:1347#0 0 chr22 31656017 37 34M * 0 0 GTGGAGAGTGCCTGCTCCAACTACGCCACCACTG a\\R``_^R\^_FVaNa_P\`\Wa^__M_\\a]^ XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
USI-EAS28:2:3:1230:1901#0 0 chr22 31656019 25 34M * 0 0 GGAGAGTGCCTGCTCCAACNACGCCACCACTGAG ]U]\``U]``a]aa`a[[[BBBBBBBBBBBBBBB XT:A:U NM:i:2 X0:i:1 X1:i:0 XM:i:2 XO:i:0 XG:i:0 MD:Z:19T12T1
USI-EAS28:3:93:186:495#0 16 chr22 31656020 37 34M * 0 0 GAGAGTGCCTGCTCCAACTACGCCACCACTGTGC BBBBBBBBBBBBBBB[Z_Y[_X_OH]__ZPMIV` XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
USI-EAS28:3:55:1400:1957#0 0 chr22 31656021 37 34M * 0 0 AGAGTGCCTGCTCCAACTACGCCACCACTGTGCA _]_`]a`bW\a`b`Q\aaaa^^[^[^^U`]^ZSP XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
HWI-EAS276:4:20:11:866#0 0 chr22 31656023 37 34M * 0 0 AGTGCCTGCTCCAACTACGCCACCACTGTGCAAG ]I]W^a_T^`aabb]bb^YW_aa`````\W[]]O XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
USI-EAS28:3:56:764:965#0 16 chr22 31656029 37 34M * 0 0 TGCTCCAACTACGCCACCACTGTGCAAGTGAAAG ][a_a_a\^W`a`baZbbaTa]`abb^SaTbbb_ XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
USI-EAS28:2:7:707:1027#0 0 chr22 31656030 37 34M * 0 0 GCTCCAACTACGCCACAACTGTGCAAGTGAAAGA ZVWY`bb]R`^T]_aVD[```^_N`aS[[\]^ZB XT:A:U NM:i:1 X0:i:1 X1:i:0 XM:i:1 XO:i:0 XG:i:0 MD:Z:16C17
USI-EAS28:2:38:109:1946#0 16 chr22 31656031 37 34M * 0 0 CTCCAACTACGCCACCACTGTGCAAGTGAAAGAG \X`_Z``^UaL`aUaa^baUaaaZ__aI_ababa XT:A:U NM:i:0 X0:i:1 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:34
Then the bed file "test.bed" is like below, only one record:
chr22 31656041 31656041
Then I run:
samtools stats --target-regions ./test.bed ./test_subset.sam
Below is the output of the above command:
# This file was produced by samtools stats (1.22.1+htslib-1.22.1) and can be plotted using plot-bamstats
# This file contains statistics for all reads.
# The command line was: stats --target-regions ./test.bed ./test_subset.sam
# CHK, Checksum [2]Read Names [3]Sequences [4]Qualities
# CHK, CRC32 of reads which passed filtering followed by addition (32bit overflow)
CHK 26fde71f 57af31e8 868ff8b1
# Summary Numbers. Use `grep ^SN | cut -f 2-` to extract this part.
SN raw total sequences: 15 # excluding supplementary and secondary reads
SN filtered sequences: 0
SN sequences: 15
SN is sorted: 1 # sorted by coordinate
SN 1st fragments: 15
SN last fragments: 0
SN reads mapped: 15
SN reads mapped and paired: 0 # paired-end technology bit set + both mates mapped
SN reads unmapped: 0
SN reads properly paired: 0 # proper-pair bit set
SN reads paired: 0 # paired-end technology bit set
SN reads duplicated: 0 # PCR or optical duplicate bit set
SN reads MQ0: 0 # mapped and MQ=0
SN reads QC failed: 0
SN non-primary alignments: 0
SN supplementary alignments: 0
SN total length: 510 # ignores clipping
SN total first fragment length: 510 # ignores clipping
SN total last fragment length: 0 # ignores clipping
SN bases mapped: 510 # ignores clipping
SN bases mapped (cigar): 168 # more accurate
SN bases trimmed: 0
SN bases duplicated: 0
SN mismatches: 6 # from NM fields
SN error rate: 3.571429e-02 # mismatches / bases mapped (cigar)
SN average length: 34
SN average first fragment length: 34
SN average last fragment length: 0
SN maximum length: 34
SN maximum first fragment length: 34
SN maximum last fragment length: 0
SN average quality: 55.6
SN insert size average: 0.0
SN insert size standard deviation: 0.0
SN inward oriented pairs: 0
SN outward oriented pairs: 0
SN pairs with other orientation: 0
SN pairs on different chromosomes: 0
SN percentage of properly paired reads (%): 0.0
SN bases inside the target: 1
SN percentage of target genome with coverage > 0 (%): 100.00
The "bases mapped (cigar)" value is 168? How this is calculated?
Below is the IGV plot, it only contains these 15 reads:
• 228 views
•
link
0 answers
No answers yet.
Log in to answer this question.