Question: How To Filter Out Sequences With No Hits From Bwa Sam Output?
3
gravatar for 2184687-1231-83-
8.1 years ago by
2184687-1231-83-5.0k wrote:

If I run bwa bwasw like this:

~/src/bwa/latest/bwa/bwa bwasw -z 1000 ~/refs/human_g1k_v37.fasta cluster-14263.walks

I get back the positive hits from some of the sequences, but also the empty hits for the sequences that don't map anywhere:

cluster-14263-walk-35   16      4       44630869        0       3S3M3D10M2I6M9D18M1I12M5I22M2I4M26S     *       0       0       CAAGTACAGCTTTTCACTGAGGCATGGAATAAATATCCATTATGGCTCTGAGTGCCCCCAGAATGCACTGCGCACCAGCGCCATGTTTTAGCTCAAGCATAACAACTCCTGACC    *       AS:i:31 XS:i:31 XF:i:1  XE:i:0  XN:i:0
cluster-14263-walk-36   4       *       0       0       *       *       0       0       GGTCAGGAGTTGTTATGCTTGAGCTAAAACATGGCGCTGGTGCGCAGTGCATTCTGGGGGCACTCAGAGCCATAATGGATATTTATTCCATGCCTCAGTGGAAG        *
cluster-14263-walk-37   4       *       0       0       *       *       0       0       GGTCAGGAGTTGTTATGCTTGAGCTAAAACATGGCGCTGGTGCGCAGTGCATTCTGGGGGCAT *

What is the cleanest/recommended way of saving only the hits and not the sequences that have no hit?

Sorry if this has been answered before, I couldn't find an answer other than grepping on the resulting files.

filter sam bwa • 2.3k views
ADD COMMENTlink modified 8.1 years ago by Docroberson280 • written 8.1 years ago by 2184687-1231-83-5.0k
6
gravatar for Aaronquinlan
8.1 years ago by
Aaronquinlan11k
United States
Aaronquinlan11k wrote:

The FLAG field has a bit for testing whether or not a read was "mapped". One can exclude unmapped reads as follows:

samtools view -F 0x4 aln.bam
ADD COMMENTlink written 8.1 years ago by Aaronquinlan11k
2

Careful, MAPQ == 0 indicates the read could not be mapped unambiguously to the ref genome.

ADD REPLYlink written 8.1 years ago by Drio910

~/bin/bwa bwasw -z 1000 $ref $reads | ~/bin/samtools view -F 0x4 -bt $ref.fai - > file.bam # Perfect! Thanks Aaron!

ADD REPLYlink written 8.1 years ago by 2184687-1231-83-5.0k

Good point drio. If you want to do "OR" logic such as this (-F 0x4 || MAPQ != 0), I would recommend bamtools' JSON "rules" files.

ADD REPLYlink written 8.1 years ago by Aaronquinlan11k
1
gravatar for Docroberson
8.1 years ago by
Docroberson280
the lab
Docroberson280 wrote:

You can also use Picard tools. http://goo.gl/ndzYU Something like: java -jar picard/directory/ViewSam.jar INPUT=aln.bam ALIGNMENT_STATUS=Unaligned

ADD COMMENTlink written 8.1 years ago by Docroberson280
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1098 users visited in the last hour