Question: HISAT2 output direct to bam
0
gravatar for Adrian Pelin
19 months ago by
Adrian Pelin2.4k
Canada
Adrian Pelin2.4k wrote:

Hello,

HISAT2 can only output sam files, which can be quite large. I found an option to output only reads that mapped to the reference but still the files can be 100s of GB.

Is there any issue in directly converting to bam as such:

${hisat2}/hisat2 -p 4 --rg-id=${4} -x $l --dta ${strand} --no-unal -1 $2 -2 $3 -U $4 | samtools view -Sbh > hisat2_output.bam

I am specifically asking about piping stdout to samtools view to convert to bam.

Thanks.

rna-seq samtools • 3.1k views
ADD COMMENTlink modified 19 months ago by ATpoint41k • written 19 months ago by Adrian Pelin2.4k
5
gravatar for ATpoint
19 months ago by
ATpoint41k
Germany
ATpoint41k wrote:

That will work fine (and is even the recommended workflow as SAM files are uncompressed and only take up space). You can even have it more efficient by doing:

hisat2 (options)... | \
  tee >(samtools flagstat - > hisat2_output.flagstat) | \
  samtools sort -O BAM | \
  tee hisat2_output.bam | \
  samtools index - hisat2_output.bam.bai

That will give you the sorted BAM, the flagstat summary statistics and the BAM index all from hisat2 stdout. tee duplicates the stream from stdin and here in combination with process substitution can be exploited to run multiple commands on essentially the same file.

ADD COMMENTlink modified 19 months ago • written 19 months ago by ATpoint41k

Ok, it makes sense to pipe to sort like swbarnes2 and yourself suggested. What is the impact on RAM usage? Piping to unsorted bam can be done on the go. Do you need to load the whole output into memory in order for you to sort?

ADD REPLYlink modified 19 months ago by ATpoint41k • written 19 months ago by Adrian Pelin2.4k
1

sort has a -m option that specifies the amount of memory to be used before it spills data as intermediate/tmp files to disk. Default is 768Mb I think per thread.

ADD REPLYlink written 19 months ago by ATpoint41k

Hi, I have a problem about this command. It's can work when I use a bash script to execute it but not in command line. When I use command line, I got the error message: -sh: syntax error near unexpected token `(' , I can't figure out what happened. Could anybody know what wrong with this?

hisat2 -q --time --novel-splicesite-outfile wt_2.ns.tsv --summary-file wt_2.sf.txt --met-file wt_2.mf.txt --threads 8 -x genome_snp_tran -1 wt_2_1P.fastq.gz -2 wt_2_2P.fastq.gz | tee >(samto
ols flagstat - > wt_2.flagstat) | samtools sort -O BAM | tee wt_2.bam | samtools index - wt_2.bam.bai
ADD REPLYlink modified 8 weeks ago by ATpoint41k • written 8 weeks ago by afukada0

Let me answer this question myself. I find that because of the "sh" issue. If you want to do this, you should be sure to use one of bash/ksh/zsh. Using for instance "sh" won't work.

ADD REPLYlink modified 8 weeks ago • written 8 weeks ago by afukada0
1
gravatar for swbarnes2
19 months ago by
swbarnes29.1k
United States
swbarnes29.1k wrote:

That should work. Better, pipe straight to samtools sort.

ADD COMMENTlink written 19 months ago by swbarnes29.1k
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: 1309 users visited in the last hour