Question: Pipe SamToFastq, BWA mem not working
0
gravatar for abbysue
8 months 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 • 616 views
ADD COMMENTlink modified 4 months ago by Marco Pessoa0 • written 8 months 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 8 months ago by Pierre Lindenbaum127k

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 8 months ago by aalith10

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 8 months ago by ATpoint31k

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 8 months ago by aalith10
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 8 months ago by WouterDeCoster43k
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 8 months ago by abbysue10

I had the same question!!

ADD REPLYlink written 8 months ago by aalith10
1

Are you on Mac?

ADD REPLYlink modified 8 months ago • written 8 months ago by ATpoint31k

Nope (this 20 character min... damn)

ADD REPLYlink written 8 months 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 8 months ago by ATpoint31k
0
gravatar for Marco Pessoa
4 months ago by
Brazil
Marco Pessoa0 wrote:

What fixed this error for me was to, in addition to indexing the reference with bwa, also indexing the reference using samtools faidx, and creating a sequence dictionary using Picard's CreateSequenceDictionary as a final step.

This (old) tutorial at gatk's forum shows details.

ADD COMMENTlink written 4 months ago by Marco Pessoa0
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: 774 users visited in the last hour