Python script for bowtie2 alignment

[mavadhani]

Hello All!

I am new to coding and am trying to write a python script to analyze a bunch of fastq files from RNASeq analysis. I am using the subprocess module to call bowtie2 in order to index the genome that is in fasta format. That works just fine. But I am stuck in the next step which is using bowtie2 to align the unpaired short reads to the genome. I know I am using the correct command as it works at the terminal generating a .sam file. But my python code is not working. Here are some of the ones I have tried--

#process = subprocess.call (['bowtie2_align', indexfiles, fastqFile1, Outputfile])
#process = subprocess.Popen(['bowtie2_align', stdin= indexfile, fastqFile1, stdout=Outputfile])

#process2 = subprocess.call(["bowtie2_align" (-x indexfile, -U fastqFile1, -S 'Outputfile.sam')],shell=True)

align = 'bowtie2_align' -x cellvibrio_index -U fastqFile1 -S initial_aligns.sam'

I would really appreciate any help/suggestions to fix this code.

Thanks,
-M

[dpryan]

It'd be something like:

subprocess.call(["bowtie2", "-x", indexBaseName, "-U", fastqFile1, "S", "Outpufile.sam"])

However, I highly encourage you to not do this and instead use something like snakemake. This lets you do things directly in python when you need to but also abstracts all of that away when you don't (like this example). Plus, you can use a cluster and in general more easily control things without reinventing the wheel for the umpteenth time.

[mavadhani]

Thank you very much dpryan. I will certainly give snakemake a shot.

-Madhavi