samtools piping with awk/ bash commands -> Wonky things happen!!!
My problem is so weird.
I'm using samtools with two windows like
$ samtools view in.bam 1:1000-2000 1:3000-4000
And it will print out paired ends with the first mate mapping to the first window and second window
When I do this
$ samtools view in.bam 1:1000-2000 1:3000-4000 | sort
or this
$ samtools view in.bam 1:1000-2000 1:3000-4000 | sort | uniq
It only prints out the first window?
When I do this
$ samtools view in.bam 1:1000-2000 1:3000-4000 | awk '{ if($7 == "=") { print $0 }}'
It will print the first window as such
HS2000-123 145 1 1005 60 100M ...
and the second ...
HS2000-456 145 1 3509 60 100M ...
The first is tab delimited and the second is space delimited
When I pipe on sort | uniq to the awk command I only return the first window.
Am I going crazy? It feels like it. I could fix this is perl, but want to know why.
What the phread is going on here?
• 3,515 views
•
link
updated
by
Ram
•
written
by
QVINTVS_FABIVS_MAXIMVS
•
written
by
QVINTVS_FABIVS_MAXIMVS
0 answers
No answers yet.
Log in to answer this question.
More posts like this
-
match bin number with the coordinate for each chromosome
written by mthm •
-
samtools view not removing unwanted headers
written by bpf24 •
-
Plotting CpG methylation
written by chaudharyc61 •
-
Creating bed files for consecutive genomic windows out from polymorphisms
written by spiral01
-
Bash- Sort a file based on a list in another file
written by jblommaert92
-
How do I replace a value used in script with a range of values in a file?
written by josephmccarter1
-
Optimizing samtools shell scripts. How to make code run faster
written by QVINTVS_FABIVS_MAXIMVS
-
What does this mean in my read (samtools)
written by QVINTVS_FABIVS_MAXIMVS
-
Plotting Fold change (FC) inside genomic interval
written by viniciushs88
-
Comparing Two Tables And Display...
written by Empyrean •
I can't reproduce this. I get both intervals regardless of piping
Do you see the same behavior if you run
\samtoolsinstead ofsamtools?If not, maybe some pre-defined aliases are disrupting your workflow.
This is very good suggestion.
Thank you :)
Thank you!
This did work, however I'm still getting this weird delimiter issue where some reads have tabs and others have spaces
So, this processes both files and writes out all windows but you still see the delimiter issue, is it? (Just clearing stuff up)
Put it into a file and process the file. Then you can test what happens - is it the samtools or the other commands that misbehave? You should not fix it in perl, this should work on its own
It is always the simplest of mistakes that produce the most mind-bendingly mysterious errors
The point of me piping commands is so I can do everything without writing to files. I have over 1.5 million regions in 250 genomes to look for mate !
The point of writing to files here is to drill down and figure out what's wrong with the pipeline, so you can correct it and get the intended use out of the pipeline. That's how testing works.