Running mutltiple Fastq files with Hisat2 in a loop

I am trying to loop paired end reads in hisat2 with an output of sorted bam using a pipe.

I am a novice at coding but am keen to improve my skills by learning how to loop.

I have created this loop using the answers from Run mutltiple Fastq files on Hisat2

 for f in `ls *.fastq.gz | sed 's/_[12].fastq.gz//g' | sort -u`
do
echo hisat2 -x /mnt/z/NathanHaffordTear/RNAseqproject/RNAseqanalysis/hg38_index -1 ${f}_1.fastq.gz -2 ${f}_2.fastq.gz | samtools sort -o ${f}.hisat.sorted.bam
done

For some reason the loop does work without the pipe but not with. How can i configure this loop to work with the piping as I do not want to have the intermediate files?

Many thanks in advance

rna-seq

• 3,820 views • link

updated 7.1 years ago by ATpoint

• written 7.1 years ago by n.tear •

Is the echo supposed to precede the call to hisat2?

• 0 views • link

written 7.1 years ago by russhh

I added it in to check the code before running.

• 0 views • link

written 7.1 years ago by n.tear •

$ ls *.gz
test_1.fastq.gz  test_2.fastq.gz

for f in `ls *.fastq.gz | sed 's/_[12].fastq.gz//g' | sort -u`
do
echo "hisat2 -x /mnt/z/NathanHaffordTear/RNAseqproject/RNAseqanalysis/hg38_index -1 ${f}_1.fastq.gz -2 ${f}_2.fastq.gz | samtools sort -o ${f}.hisat.sorted.bam"
done

output:

hisat2 -x /mnt/z/NathanHaffordTear/RNAseqproject/RNAseqanalysis/hg38_index -1 test_1.fastq.gz -2 test_2.fastq.gz | samtools sort -o test.hisat.sorted.bam

OP is supposed to run:

for f in `ls *.fastq.gz | sed 's/_[12].fastq.gz//g' | sort -u`
do
echo "parising f...wait..." && hisat2 -x /mnt/z/NathanHaffordTear/RNAseqproject/RNAseqanalysis/hg38_index -1 ${f}_1.fastq.gz -2 ${f}_2.fastq.gz | samtools sort -o ${f}.hisat.sorted.bam 
done

OP can also try parallel:

parallel --plus --dry-run  "hisat2 -x /mnt/z/NathanHaffordTear/RNAseqproject/RNAseqanalysis/hg38_index -1 {} -2 {=s/1/2/=} | samtools sort -o {=s/_1.fastq.gz//=}.hisat.sorted.bam" ::: *_1.fastq.gz

output:

hisat2 -x /mnt/z/NathanHaffordTear/RNAseqproject/RNAseqanalysis/hg38_index -1 test_1.fastq.gz -2 test_2.fastq.gz | samtools sort -o test.hisat.sorted.bam

• 0 views • link

written 7.1 years ago by cpad0112

1 answer

Remove the echo as you want to execute rather than print the command. If you want "status updates" do

echo $f && hisat2 (...)

• 0 views • link

written 7.1 years ago by ATpoint

Hi thanks for the reply

Im aware of this but want to make sure the code works with echo before running

running my code I get the error:

[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting
[W::sam_read1] Parse error at line 1
samtools sort: truncated file. Aborting

I will try with the added echo $f && hisat2 (...) and let you know how it goes

thanks

• 0 views • link

written 7.1 years ago by n.tear •

Of course because hisat2 does not produce output and therefore samtools does not get any valid stdin. Remove the echo and maybe put a head like hisat2 | head -n 1000 | samtools sort to get an idea if it works properly (and it should, the command is ok).

• 0 views • link

written 7.1 years ago by ATpoint

Thanks that makes perfect sense. doh

The command you suggested before seems to be working great, i like getting the updates

• 0 views • link

written 7.1 years ago by n.tear •