Question: How to get unique mapped reads?
0
gravatar for biolab
3.7 years ago by
biolab1.1k
biolab1.1k wrote:

Dear all,

Previously I thought when the -g/--max-multihits option in tophat was set to 1, then I acquired the unique mapped reads.  However, after careful thinking, I am probably wrong. 

Now my idea is: when doing mapping, set the -g/--max-multihits option to 2 or above 2, and then write a script to eliminate those reads occurring twice or more times.  Am I right?  Are there any tools or simple approach to obtain unique mapped reads?   I appreciate any of your comments and answers.  THANKS A LOT!

mapping bam • 3.1k views
ADD COMMENTlink modified 3.7 years ago by QVINTVS_FABIVS_MAXIMVS2.2k • written 3.7 years ago by biolab1.1k
1

Using BBMap, you can set the flag "ambig=toss" which will ensure all mapped reads have unique mappings.  Reads that multimap will instead be marked as unmapped (and discarded, if you use "outm" instead of "out").

ADD REPLYlink written 3.7 years ago by Brian Bushnell16k
1

You can fgrep all the reads with tag NH:i:1, which filter outs only unique mapped reads.

sam/bam | fgrep -w NH:i:1

ADD REPLYlink written 3.7 years ago by tharveshliyakat60
2
gravatar for Philipp Bayer
3.7 years ago by
Philipp Bayer6.0k
Australia/Perth/UWA
Philipp Bayer6.0k wrote:

This has been asked before a few times - here's one discussion about bowtie2 (since tophat uses bowtie2 internally you could use the same points)

How to extract unique mapped results from Bowtie2 bam results?

Most notably:

>You're better served by simply filtering on a meaningful MAPQ (5 or 10 are often reasonable choices), which has the benefit of actually doing what you want, namely filtering according the the likelihood that an alignment is correct.

So using "samtools view -q 10 input > filtered_output" (or even higher? check your distribution of scores)

ADD COMMENTlink written 3.7 years ago by Philipp Bayer6.0k
2
gravatar for QVINTVS_FABIVS_MAXIMVS
3.7 years ago by
USA SoCal
QVINTVS_FABIVS_MAXIMVS2.2k wrote:
    \samtools view -f 0x02 in.bam | sort | uniq >unique_mapped_reads.txt

Obviously you would want to specify a region

    \samtools view -f 0x02 in.bam 15:20000000-30000000 | sort | uniq >unique_mapped_praderWilliLocus_reads.txt

for example.

ADD COMMENTlink modified 3.7 years ago • written 3.7 years ago by QVINTVS_FABIVS_MAXIMVS2.2k
1

Thank you all for your comments and answers.  They are definitely helpful for my work.

ADD REPLYlink written 3.7 years ago by biolab1.1k

you work on CNVs too mate ;)

ADD REPLYlink written 3.7 years ago by QVINTVS_FABIVS_MAXIMVS2.2k
1
gravatar for seidel
3.7 years ago by
seidel6.8k
United States
seidel6.8k wrote:

On thing you could do is to examine the NH tag in the tophat output file (accepted_hits.bam). It lists the number of alignments for that read, and thus you could filter out reads with NH > 1.

ADD COMMENTlink written 3.7 years ago by seidel6.8k
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: 1144 users visited in the last hour