<p>Other options... <a href="http://www.htslib.org/doc/samtools.html" rel="nofollow noreferrer">samtools stats</a> maybe gives you some of the stats you are after. Chances are you already have samtools so nothing additional to install. Also useful may be <a href="https://github.com/alimanfoo/pysamstats" rel="nofollow noreferrer">pysamstats</a>, in particular the <code>variation</code> sub-command.</p>

CC BY-SA 4.0
2018-11-21T15:50:55.440Z
pysamstats looks really neat!
winni2k
CC BY-SA 3.0
2018-02-03T15:37:25.633Z
2018-02-03T15:37:25.633Z
dariober
<p>I've written my own script, <a href="https://github.com/gringer/bioinfscripts/blob/master/readstomper.pl" rel="nofollow noreferrer">readstomper</a>, which gives per-site statistics from mpileup results. Here's an example usage:</p>

<pre><code>$ samtools mpileup -B -Q 0 -f circ-Nb-ec3-mtDNA.fasta LAST_OL_132394_vs_mtDNA.bam | \
    /bioinf/scripts/readstomper.pl &gt; LAST_OLmpileupProp_132394_vs_mtDNA.csv

Assembly,Position,Coverage,ref,cR,pR,A,C,G,T,d,i,InsMode
Nb_mtDNA,11456,86,T,79.00,0.92,0.0000,0.0465,0.0000,0.0000,0.03,0.01,A;1.00
Nb_mtDNA,11457,86,A,81.00,0.94,0.0000,0.0000,0.0233,0.0000,0.03,0.02,TA;0.50
Nb_mtDNA,11458,86,G,79.00,0.92,0.0349,0.0233,0.0000,0.0116,0.01,0.00
Nb_mtDNA,11459,86,A,78.00,0.92,0.0000,0.0000,0.0235,0.0471,0.01,0.00
Nb_mtDNA,11460,86,T,5.00,0.06,0.0233,0.8140,0.0233,0.0000,0.08,0.00
Nb_mtDNA,11461,86,G,57.00,0.66,0.1512,0.1163,0.0000,0.0465,0.02,0.01,AATAACAACACGTAACCGA;1.00
Nb_mtDNA,11462,86,G,78.00,0.91,0.0349,0.0349,0.0000,0.0116,0.01,0.00
Nb_mtDNA,11463,86,G,66.00,0.77,0.0465,0.0698,0.0000,0.0116,0.10,0.01,GCAATA;1.00
Nb_mtDNA,11464,86,T,74.00,0.86,0.0000,0.0814,0.0233,0.0000,0.03,0.00
</code></pre>

<p>It defaults to reporting variants as proportions of the total coverage, but this can be changed to reporting the individual counts:</p>

<pre><code>$ samtools mpileup -B -Q 0 -f mmus_chrM.fa mm2_1.2M_vs_Mmus_mtDNA.bam | \
    /bioinf/scripts/readstomper.pl --counts

Assembly,Position,Coverage,ref,skip,cR,pR,A,C,G,T,d,i,InsMode
chrM,1528,491,C,0,409,0.87,17,0,2,15,25,65,A;52
chrM,1529,492,A,0,366,0.75,0,2,2,9,112,12,A;5
chrM,1530,493,A,0,409,0.83,0,2,1,13,68,7
chrM,1531,495,A,0,463,0.94,0,1,10,0,21,8,G;5
chrM,1532,496,A,0,463,0.94,0,2,5,8,16,4,G;3
chrM,1533,495,G,0,459,0.93,3,1,0,2,30,4
chrM,1534,496,A,0,451,0.91,0,0,4,1,39,23,G;21
chrM,1535,503,G,0,391,0.78,4,1,0,4,103,9,A;6
chrM,1536,506,G,0,496,0.98,6,0,0,0,4,10,A;3
chrM,1537,506,G,0,488,0.96,8,1,0,4,5,5,A;3
chrM,1538,506,A,0,483,0.95,0,4,3,4,12,12
chrM,1539,506,C,0,466,0.92,9,0,2,14,14,27,A;6
chrM,1540,506,A,0,469,0.93,0,0,17,6,14,11,T;3
chrM,1541,506,G,0,453,0.90,12,1,0,5,31,21,A;4
chrM,1542,505,C,0,421,0.83,6,0,2,14,62,8,CT;3
</code></pre>

CC BY-SA 3.0
2018-02-02T21:29:17.937Z
2018-02-02T21:29:17.937Z
gringer
<p>using BioAlcidaeJdk: <a href="http://lindenb.github.io/jvarkit/BioAlcidaeJdk.html" rel="nofollow noreferrer">http://lindenb.github.io/jvarkit/BioAlcidaeJdk.html</a> , you can loop over the cigar strings of each mapped read and extract the number of inserted bases, clipped bases, etc...</p>

<pre><code>java -jar dist/bioalcidaejdk.jar -e 'stream().filter(R-&gt;!R.getReadUnmappedFlag() &amp;&amp; R.getCigar()!=null).flatMap(R-&gt;R.getCigar().getCigarElements().stream()).collect(Collectors.toMap(C-&gt;C.getOperator(),C-&gt;(long)C.getLength(),(L1,L2)-&gt;L1+L2)).forEach((O,N)-&gt;println(O.name()+"\t"+N));' toy.bam
I   15
H   11
S   1
D   1
N   14
P   2
M   191
</code></pre>

<p>to get the number of mismatches, you can preprocess the bam with <a href="http://lindenb.github.io/jvarkit/SamFixCigar.html" rel="nofollow noreferrer">http://lindenb.github.io/jvarkit/SamFixCigar.html</a> </p>

CC BY-SA 3.0
2018-02-02T21:34:06.243Z
2018-02-02T21:34:06.243Z
Pierre
<p>I used BWA-MEM to alignment and I would like to gather the some informations like total % of match, mismatch, insert/delete etc. I am wondering if there is any existing tools that produces this information?</p>

<phylogenetics><samtools><sam>
Parsing SAM/BAM files for additional information
# Parsing SAM/BAM files for additional information

(Asked by: username_0 on 2018-02-02)

I used BWA-MEM to alignment and I would like to gather the some informations like total % of match, mismatch, insert/delete etc. I am wondering if there is any existing tools that produces this information?

**Tags:** phylogenetics, samtools, sam

---

## 3 answers

### Answer by username_1 on 2018-02-03. Score: 2

Other options... [samtools stats](http://www.htslib.org/doc/samtools.html) maybe gives you some of the stats you are after. Chances are you already have samtools so nothing additional to install. Also useful may be [pysamstats](https://github.com/alimanfoo/pysamstats), in particular the `variation` sub-command.

> **Comments:**
> **username_2**: pysamstats looks really neat!  

---

### Answer by username_3 on 2018-02-02. Score: 0

I've written my own script, readstomper, which gives per-site statistics from mpileup results. Here's an example usage:
$ samtools mpileup -B -Q 0 -f circ-Nb-ec3-mtDNA.fasta LAST_OL_132394_vs_mtDNA.bam | \
    /bioinf/scripts/readstomper.pl > LAST_OLmpileupProp_132394_vs_mtDNA.csv

Assembly,Position,Coverage,ref,cR,pR,A,C,G,T,d,i,InsMode
Nb_mtDNA,11456,86,T,79.00,0.92,0.0000,0.0465,0.0000,0.0000,0.03,0.01,A;1.00
Nb_mtDNA,11457,86,A,81.00,0.94,0.0000,0.0000,0.0233,0.0000,0.03,0.02,TA;0.50
Nb_mtDNA,11458,86,G,79.00,0.92,0.0349,0.0233,0.0000,0.0116,0.01,0.00
Nb_mtDNA,11459,86,A,78.00,0.92,0.0000,0.0000,0.0235,0.0471,0.01,0.00
Nb_mtDNA,11460,86,T,5.00,0.06,0.0233,0.8140,0.0233,0.0000,0.08,0.00
Nb_mtDNA,11461,86,G,57.00,0.66,0.1512,0.1163,0.0000,0.0465,0.02,0.01,AATAACAACACGTAACCGA;1.00
Nb_mtDNA,11462,86,G,78.00,0.91,0.0349,0.0349,0.0000,0.0116,0.01,0.00
Nb_mtDNA,11463,86,G,66.00,0.77,0.0465,0.0698,0.0000,0.0116,0.10,0.01,GCAATA;1.00
Nb_mtDNA,11464,86,T,74.00,0.86,0.0000,0.0814,0.0233,0.0000,0.03,0.00

It defaults to reporting variants as proportions of the total coverage, but this can be changed to reporting the individual counts:
$ samtools mpileup -B -Q 0 -f mmus_chrM.fa mm2_1.2M_vs_Mmus_mtDNA.bam | \
    /bioinf/scripts/readstomper.pl --counts

Assembly,Position,Coverage,ref,skip,cR,pR,A,C,G,T,d,i,InsMode
chrM,1528,491,C,0,409,0.87,17,0,2,15,25,65,A;52
chrM,1529,492,A,0,366,0.75,0,2,2,9,112,12,A;5
chrM,1530,493,A,0,409,0.83,0,2,1,13,68,7
chrM,1531,495,A,0,463,0.94,0,1,10,0,21,8,G;5
chrM,1532,496,A,0,463,0.94,0,2,5,8,16,4,G;3
chrM,1533,495,G,0,459,0.93,3,1,0,2,30,4
chrM,1534,496,A,0,451,0.91,0,0,4,1,39,23,G;21
chrM,1535,503,G,0,391,0.78,4,1,0,4,103,9,A;6
chrM,1536,506,G,0,496,0.98,6,0,0,0,4,10,A;3
chrM,1537,506,G,0,488,0.96,8,1,0,4,5,5,A;3
chrM,1538,506,A,0,483,0.95,0,4,3,4,12,12
chrM,1539,506,C,0,466,0.92,9,0,2,14,14,27,A;6
chrM,1540,506,A,0,469,0.93,0,0,17,6,14,11,T;3
chrM,1541,506,G,0,453,0.90,12,1,0,5,31,21,A;4
chrM,1542,505,C,0,421,0.83,6,0,2,14,62,8,CT;3


---

### Answer by username_4 on 2018-02-02. Score: 0

using BioAlcidaeJdk: http://lindenb.github.io/jvarkit/BioAlcidaeJdk.html , you can loop over the cigar strings of each mapped read and extract the number of inserted bases, clipped bases, etc...
java -jar dist/bioalcidaejdk.jar -e 'stream().filter(R->!R.getReadUnmappedFlag() && R.getCigar()!=null).flatMap(R->R.getCigar().getCigarElements().stream()).collect(Collectors.toMap(C->C.getOperator(),C->(long)C.getLength(),(L1,L2)->L1+L2)).forEach((O,N)->println(O.name()+"\t"+N));' toy.bam
I   15
H   11
S   1
D   1
N   14
P   2
M   191

to get the number of mismatches, you can preprocess the bam with http://lindenb.github.io/jvarkit/SamFixCigar.html


---

