Entering edit mode
8.6 years ago
ddzhangzz
▴
90
I was trying to align paired end RNASeq data to a locally-built reference using bowtie2:
$bowtie2 --trim3 0 --threads 6 --no-discordant --no-mixed -k 1 \
-x $idx -1 $f1 -2 $f2 -S classI.sam
And the sam file output looks like this (ignored the header):
ERR356371.7445 141 * 0 0 * * 0 0 ATTGAGAAGCCCCCCTTCGAGCTGCCAGACTTCATCAAACGCACAGGCA 2@@@DDDDA?DFDHGBHHICEGGDFDGCCGIIIJJBFHD<DFHIJAAGH YT:Z:UP
ERR356371.7446 77 * 0 0 * * 0 0 CTACCATCACCCTCTACATCACCGCCCCGACCTTAGCTCTCACCATCGC CCCFFFFFHHHHHJJIJJJJIJJIJJJJJIJJJJJJJJIJGGIJEHJJI YT:Z:UP
ERR356371.7446 141 * 0 0 * * 0 0 ATCGCAGTGCGCCGATCAGGGCGTAGTTTGAGTTTGATGCTCACCCTGA =CCCFFFDFHHHHHGIJJJJJJJHIGFIIJGIIIJJIJJJJJJJJJJJH YT:Z:UP
.
.
ERR356371.7400 99 HLA:HLA05978 527 255 49M = 636 158 CGGCGGAGCAGCGGAGAGTCTACCTGGAGGGCCGGTGCGTGGACGGGCT CCCFFFFFGHHGHJJIIJHFHJIIJJJJJJJJJHF=BDDBDDDDDDDDD AS:i:0 XN:i:0 XM:i:0 XO:i:0 XG:i:0 NM:i:0 MD:Z:49 YS:i:-8 YT:Z:CP
ERR356371.7400 147 HLA:HLA05978 636 255 49M = 527 -158 TATGACCCACCACCCCATCTCTGACCATGAGGCCACCCTGAGGTGCAGT JJIHGJHF@HD?JJIJIHHGJIIGJJJIJJJJIGEHHHHHFFFFFCCC2 AS:i:-8 XN:i:0 XM:i:2 XO:i:0 XG:i:0 NM:i:2 MD:Z:46T1G0 YS:i:0 YT:Z:CP
.
.
.
ERR356371.7447 77 * 0 0 * * 0 0 AATTCATGTCTGCACCAACTTGGACAACTGGGGGAGGCTGGATGCATCA CCCFFFFFHHHHHJJJJJJJJJJIJJJJJJJJJJGHJIJIJGJJJIJJI YT:Z:UP
ERR356371.7447 141 * 0 0 * * 0 0 ATTGAAGTCAAGGTGTTGGCTGCCACATGCCATGTCCTCTTAACAAGAT =CCCFFFFFHHHHFGIJJJJJJJJIJJJJJJIJJIJJJJJJJJJJJJJJ YT:Z:UP
ERR356371.7448 77 * 0 0 * * 0 0 AATTGCCAATATTTGCCCCTCCCCTTTAATAAAACTTTTAAGAAGTCAC <?@DDFFDHFFFHJGEFBFFGAEDGHGIJIFIEF@GHGHBHEBFHGFIH YT:Z:UP
ERR356371.7448 141 * 0 0 * * 0 0 AGGCAAATATTGGCAATTAGTTGGCCGTAAGTACAGCACAGTGCAGCTT A??@D?DDDDFHDH@GGGIGAACDBHGCCEFIIIIG@DDDBHGGGGGHI YT:Z:UP
My question is whether there is an option in bowtie2 only output the valid alignments such as:
ERR356371.7400 99 HLA:HLA05978 527 255 49M = 636 158 CGGCGGAGCAGCGGAGAGTCTACCTGGAGGGCCGGTGCGTGGACGGGCT CCCFFFFFGHHGHJJIIJHFHJIIJJJJJJJJJHF=BDDBDDDDDDDDD AS:i:0 XN:i:0 XM:i:0 XO:i:0 XG:i:0 NM:i:0 MD:Z:49 YS:i:-8 YT:Z:CP
ERR356371.7400 147 HLA:HLA05978 636 255 49M = 527 -158 TATGACCCACCACCCCATCTCTGACCATGAGGCCACCCTGAGGTGCAGT JJIHGJHF@HD?JJIJIHHGJIIGJJJIJJJJIGEHHHHHFFFFFCCC2 AS:i:-8 XN:i:0 XM:i:2 XO:i:0 XG:i:0 NM:i:2 MD:Z:46T1G0 YS:i:0 YT:Z:CP
Namely the middle part of the full results. Does somebody know how to output this?
it does make more sense not to print out unaligned reads if you already know that you're not interested in them, plus it would be faster. in case you aren't completely sure about it, I would suggest to convert the sam file into a bam file (
samtools sort
followed bysamtools view -Sb
) and then extract from that bam file only the aligned reads (samtools view -F4
). you will end up having to deal with temporal files though.