contig_sigs.py

import screed
import os
from os import path
import subprocess
from subprocess import Popen, PIPE

def make_dir(dirname):
    if os.path.isdir(dirname) == False:
        os.mkdir(dirname)
        print("Directory created:", dirname)

def execute(basedir,assemblies):
     for fasta_file in assemblies:
             if fasta_file.endswith('fa'):
                 genus_species = fasta_file.split('/')[5].split(".")[0]
             else:
                 genus_species = 'F_het_ref'
             print(genus_species)
             dirname = genus_species + "/"
             make_dir(dirname)
             for n, record in enumerate(screed.open(fasta_file)):
                 tmp_filename = dirname + genus_species + "." + record.name.split(" ")[0] + ".fasta"
                 #print(tmp_filename)
                 print("Computing sig for contig: "+str(n))
                 sig_filename = dirname + genus_species + "." + record.name.split(" ")[0] + ".sig"
                 #print(sig_filename)
                 with open(tmp_filename,'w') as fp:
                     name = """{}-{}""".format(genus_species,record.name)
                     fp.write(">{name}\n{seq}\n".format(name=name, seq=record.sequence))
                 sourmash_command="sourmash compute --scaled 500 "+tmp_filename+" -o "+sig_filename+" -k 51"
                 print(sourmash_command)
                 t=subprocess.Popen(sourmash_command, shell=True)
                 t.wait()
                 os.remove(tmp_filename)

basedir = "/mnt/home/ljcohen/sourmash_ont/"
A_xen = '/mnt/scratch/ljcohen/A_xenica_PromethION/A_xenica.fa'
F_not = '/mnt/scratch/ljcohen/F_nottii_PromethION/F_nottii.fa'
F_oli = '/mnt/scratch/ljcohen/F_olivaceus/F_olivaceus_masurca.fa'
F_het_ref = '/mnt/home/ljcohen/reference/GCF_000826765.1_Fundulus_heteroclitus-3.0.2_genomic.fna'
#assemblies = [A_xen,F_not,F_oli,F_het_ref]
assemblies = [F_het_ref]
execute(basedir,assemblies)