Running a hisat2 loop, samtools view and sort using python

Hi all,

I am running a large RNA seq dataset and have the trimmed reads ready for alignment through hisat2. If I try and run more than one alignment at a time, then my computer will not really operate because of the strain from the function. I am newer to python and am not quite sure on how to access functions like that through it, but here is what I have so far, yet its not quite working.

import os
import glob 
for file in folder:
     if '*forward.fastq.gz' == 'reverse.fastq.gz':
       file1 ='*forward.fastq.gz'
       file2 = '*reverse.fastq.gz'
       out = file1.replace('-forward.fastq.gz', '_sorted.bam')
       os.popen('/path/to/hisat2/hisat2 -p 20  -t --fr --dta-cufflinks -x hg38_index -1 %s -2 %s | /path/to/samtools view -b | /path/to/samtools sort -o %s' %(file1, file2, out))

Thanks so much

rna-seq rna-seq python hisat2 samtools

• 432 views • link

updated 5.9 years ago by ATpoint

• written 5.9 years ago by noahhelton98 •

1 answer

Don't make it too complicated. For such a simple task a plain bash script is enough:

Say your files are called like sample1_1.fastq.gz, sample1_2.fastq.gz, sample2_1.fastq.gz, sample2_2.fastq.gz (...):

for i in *_1.fastq.gz
  do
  Basename=${i%_1.fastq.gz}
  hisat2 (options...) -x hg38_index -1 ${Basename}_1.fastq.gz -2 ${Basename}_2.fastq.gz \
  | samtools sort -o ${Basename}_sorted.bam
  done

No need for view, as sort can read SAM files directly from hisat2.

• 0 views • link

written 5.9 years ago by ATpoint

Thanks so much, I am trying to run this code, but for some reason it keeps saying their is no such directory as _1.fastq.gz and "gzip: can't stat: _2.fastq.gz (_2.fastq.gz.gz)". For my loop I had put

for i in trimmed_reads/*_1.fastq.gz 
do ... what you said had in the code

my file names are named like: iPS-CTL-Chr_S1_1.fastq.gz iPS-CTL-Chr_S1_2.fastq.gz iPS-CTL-Nuc_S2_1.fastq.gz iPS-CTL-Nuc_S2_2.fastq.gz iPS-U4ATAC-Chr_S1_1.fastq.gz iPS-U4ATAC-Chr_S1_2.fastq.gz

etc.

Any suggestions on how to fix this error? Thank you so much for your help though, I am just starting to get used to coding and bioinformatics recently, so sorry if my question is very simple.

• 0 views • link

written 5.9 years ago by noahhelton98 •