hello every one , I appear here with a serious question regarding unique mapping with bowtie2. uniqueness dos not mean to force any read to map at single place on reference chromosome [by using -k 1 in bowtie2] I have mapped my illumina reads at reference by using parameter -k 2 and -k 1
At -k 2 :
102886046 reads; of these:
102886046 (100.00%) were unpaired; of these:
29263187 (28.44%) aligned 0 times
55532464 (53.97%) aligned exactly 1 time
18090395 (17.58%) aligned >1 times
71.56% overall alignment rate
73.17% overall alignment rate
At -k 1:
102886046 (100.00%) were unpaired; of these:
29263187 (28.44%) aligned 0 times
73622859 (%) aligned exactly 1 time 0 (0%) aligned >1 times 71.56% overall alignment rate
The -k 2 parameter showing just additive results which is quite frustrating, how we would proceed if we are looking for "uniquely mapped reads".
what parameter is fully define uniquely map??????
thank you
If you want only the uniquely mapped reads from bowtie2 sam files then you can proceed as follows:
grep -E "@|NM:" bt2output.sam | grep -v "XS:" > uniq_bt2output.sam
This first checks to see if there is alignment and then removes reads with secondary alignments. Hope this helps.
If, at the end of the day, all you want are unique alignments, then don't feed bowtie2 the "-k" parameter. With default options, it will return return the best found alignment and set the XS and AS auxiliary field. For the most part, one can simply look at the MAPQ field to judge if it's a multimapper, but you can also compare the AS/XS tags (the XS tags can be missing if there are no second best alignments).
if mapping quality in 6th column of Sam file is vary from 0 to more than 42 and AS:i:0 and XS:i: has minus values like XS:i:-6 -15 and -20 etc. then how do we customize if I need MQUAL value more than 30 say for example and do we relate this MQUAL value with AS:i and XS:i
thank you
I'm not entirely sure how to parse what you wrote, so I'll rewrite what I think you're asking: "After alignment, how do we filter the results to leave only those alignments most likely to be unique?" You have another question regarding how bowtie2 calculates MAPQ, which can only really be answered by looking at the source code (a brief scroll through the source code suggests that bowtie2 has at least 3 versions of a Mapq calculator!). For filtering by MAPQ, just use "samtools view -q ...", which is the most convenient. I'm not aware of a premade script solely to filter by the AS and XS tags, but you could easily write one, though that's probably overkill.
Log in to answer this question.
