Hi
I have tried to map this short sequence
>sample_1_x12179
TCCTGTACTGAGCTGCCCCGA
which exists on chr8 twice, using various aligners.
BLAT, BLAST, and Bowtie2 found it. BWA MEM and SOAP2 did not. How come? BWA MEM finds it if I include some extra flanking bases. Is it because the sequence is too short?
BWA MEM command:
bwa mem -t 8 $bwaIndex tmp.fasta
SOAP command:
./soap2.20 \
-a tmp.fasta \
-o ${goutputfile}.tmp \
-D /mnt/NGS01/ReferenceDB/mirTools/db/genome/hsa.fa.index \
-M 0 -r 2 -v 5 -p 8 -l 15
I also tried -M 4 (find best hits), but still no results. If I remove the last two bases from the above sequence, then SOAP2 finds a hit on chr11 that requires 2 mismatches.
Thank you.
•
written
by
francisckantor •
BWA works for me (but not bwa mem, because it requires 70bp minimum query sequence length):
bwa aln $GRCh37 test.fa > test.sai
bwa samse $GRCh37 test.sai test.fa > test.sam
sample_1_x12179 16 8 41518003 0 21M * 0 0TCGGGGCAGCTCAGTACAGGA * XT:A:R NM:i:0 X0:i:2 X1:i:0 XM:i:0 XO:i:0 XG:i:0 MD:Z:21 XA:Z:8,+41517962,21M,0;
_____________________________^_______________________________^^^
The documentation for SOAP2 sucks, but I'm guessing it's an issue with your query read length.
Thank you. It seems to work. I'll read more up on the various parameter settings.
Log in to answer this question.
You're not using the
bwa memcommand in your example. And it wouldn't work anyway becausebwa memexpects a 70bp minimum query sequence length for aligning. Does runningbwa alnwithout thememlike you've shown work okay?Sorry, I did use 'bwa mem', just accidentally posted the wrong command. I'll try 'bwa aln' too. Thank you for your answer.
Can you post your SOAP command as well?
I've edited my post and added the SOAP2 command.