Question: Pipe SamToFastq, BWA mem not working
0
gravatar for abbysue
5 weeks ago by
abbysue10
abbysue10 wrote:

I'm going through the GATK pipeline to realign bam files. I'm using this command line:

gatk SamToFastq -I markedadapters.bam -F /dev/stdout --CLIPPING_ATTRIBUTE XT --CLIPPING_ACTION 2 --INTERLEAVE true -NON_PF true -TMP_DIR= path/hello | bwa mem -M -t 12 -p /mydata/Homo_sapiens_assembly38.fasta /dev/stdin | gatk MergeBamAlignment --ALIGNED_BAM /dev/stdin --UNMAPPED_BAM unmap.bam -O piped.bam -R /mydata/Homo_sapiens_assembly38.fasta --CREATE_INDEX true --ADD_MATE_CIGAR true --CLIP_ADAPTERS false --CLIP_OVERLAPPING_READS true --INCLUDE_SECONDARY_ALIGNMENTS true --MAX_INSERTIONS_OR_DELETIONS -1 --PRIMARY_ALIGNMENT_STRATEGY MostDistant --ATTRIBUTES_TO_RETAIN XS -TMP_DIR= path/hello

I keep getting this error:

htsjdk.samtools.SAMException: Error in writing fastq file /dev/stdout

I'm about to give up and just run each command separately, but these files are ~150G :(

pipe bwa sam gatk • 214 views
ADD COMMENTlink modified 29 days ago by Biostar ♦♦ 20 • written 5 weeks ago by abbysue10
3

I'm about to give up and just run each command separately,

give samtools collate | samtools fastq a try.

ADD REPLYlink written 5 weeks ago by Pierre Lindenbaum122k

I've read that the behavior of Picard in regressing to a fastq file is preferable (forum post on Biostars, have to find it real quick)

ADD REPLYlink written 5 weeks ago by aalith0

I do not see why this would be the case. I am using the above command (and never used Picard for it) all the time, never fails :)

ADD REPLYlink written 5 weeks ago by ATpoint21k

Could abbysue's problem be due to the fact that it's running inside docker and not on the local machine? the piped command listed should work

ADD REPLYlink written 5 weeks ago by aalith0
2

Could you try running the first part before the pipe and see if that works?

gatk SamToFastq -I markedadapters.bam -F /dev/stdout --CLIPPING_ATTRIBUTE XT --CLIPPING_ACTION 2 --INTERLEAVE true -NON_PF true -TMP_DIR= path/hello > mytest.fastq

If one of the steps later fails then the pipe will break. Make sure all components work before making a large piped command.

ADD REPLYlink written 5 weeks ago by WouterDeCoster40k
1

I ended up running

gatk SamToFastq -I markedadapters.bam -F f797.fastq --CLIPPING_ATTRIBUTE XT --CLIPPING_ACTION 2 --INTERLEAVE true -NON_PF true -TMP_DIR= path/hello

which of course worked. I tried your command and it worked!! So should I include the " > file.fastq" blurb to get my piped commands to work?

ADD REPLYlink written 5 weeks ago by abbysue10

I had the same question!!

ADD REPLYlink written 5 weeks ago by aalith0
1

Are you on Mac?

ADD REPLYlink modified 5 weeks ago • written 5 weeks ago by ATpoint21k

Nope (this 20 character min... damn)

ADD REPLYlink written 5 weeks ago by abbysue10
1

Ok, on Mac I sometimes had issues with non-writable stdout. As for the character limit, just add whitespaces ;-)

ADD REPLYlink written 5 weeks ago by ATpoint21k
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: 1443 users visited in the last hour