# HG changeset patch # User tyty # Date 1413783233 14400 # Node ID b0243634fe6b27c03a175326f3ecc0847712e039 # Parent 74ae6f6a3383c4fcc52f962c687d3a0a91478bea Deleted selected files diff -r 74ae6f6a3383 -r b0243634fe6b Iterative_mapping/.DS_Store Binary file Iterative_mapping/.DS_Store has changed diff -r 74ae6f6a3383 -r b0243634fe6b Iterative_mapping/iterative_map.py --- a/Iterative_mapping/iterative_map.py Mon Oct 20 01:33:11 2014 -0400 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,103 +0,0 @@ -#!/usr/bin/env python -# -*- coding: utf-8 -*- - -import sys -import os -from read_file import * -from read_s_file import * - -type_input = sys.argv[1] -seq_file = sys.argv[2] -ref_file = sys.argv[3] -shift = sys.argv[4] -length = sys.argv[5] -map_type = sys.argv[6] -output_file = sys.argv[7] - - -if map_type!="default": - s = "" - s = s+"-v "+sys.argv[8] - s = s+" -5 "+sys.argv[9] - s = s+" -3 "+sys.argv[10] - s = s+" -k "+sys.argv[11] - if sys.argv[12]: - s = s+" -a" - if int(sys.argv[13])>=1: - s = s+" -m "+sys.argv[13] - if sys.argv[14]: - s = s+" --best --strata " - -else: - s = "-v 3 -a --best --strata " - -ospath = os.path.realpath(sys.argv[0]) -ost = ospath.split('/') -syspath = "" -for i in range(len(ost)-1): - syspath = syspath+ost[i].strip() - syspath = syspath+'/' - -os.system("bowtie-build -f "+ref_file+" "+syspath+"ref > "+syspath+"log.txt") - -os.system("cp "+seq_file+" "+syspath+"seq0.fa") - -if type_input == "fasta": - tp = 'fasta' -if type_input == "fastq": - tp = 'fastq' - -k = 0 -print(type_input) -while(True): - if type_input == "fasta": - os.system("bowtie "+s+"-f "+syspath+"ref"+" "+syspath+"seq"+str(k)+".fa --quiet -S > "+syspath+"map"+str(k)+".sam") - if type_input == "fastq": - os.system("bowtie "+s+"-q "+syspath+"ref"+" "+syspath+"seq"+str(k)+".fa --quiet -S > "+syspath+"map"+str(k)+".sam") - os.system("samtools view -Sb -F 0xfff "+syspath+"map"+str(k)+".sam > "+syspath+"mapped"+str(k)+".bam 2>"+syspath+"log.txt") #get mapped reads - os.system("samtools view -Sb -f 0x4 "+syspath+"map"+str(k)+".sam > "+syspath+"umapped"+str(k)+".bam 2>"+syspath+"log.txt") #get unmapped reads - os.system("samtools view -Sb -f 0x10 "+syspath+"map"+str(k)+".sam > "+syspath+"rmapped"+str(k)+".bam 2>"+syspath+"log.txt") #get reversed mapped reads - os.system("samtools merge -f "+syspath+"unmapped"+str(k)+".bam "+syspath+"umapped"+str(k)+".bam "+syspath+"rmapped"+str(k)+".bam") #get reversed mapped reads - os.system("samtools view -h -o "+syspath+"unmapped"+str(k)+".sam "+syspath+"unmapped"+str(k)+".bam") #get reversed mapped reads - if k>0: - os.system("samtools view -h -o "+syspath+"mapped"+str(k)+".sam "+syspath+"mapped"+str(k)+".bam") #get reversed mapped reads - os.system("cut -f 1 "+syspath+"unmapped"+str(k)+".sam > "+syspath+"unmapped"+str(k)+".txt") - os.system("cut -f 1 "+syspath+"mapped"+str(k)+".sam > "+syspath+"mapped"+str(k)+".txt") - os.system("python "+syspath+"remove_map.py "+syspath+"unmapped"+str(k)+".txt "+syspath+"mapped"+str(k)+".txt "+syspath+"runmapped"+str(k)+".txt") - os.system("rm "+syspath+"mapped"+str(k)+".sam") - os.system("rm "+syspath+"mapped"+str(k)+".txt") - os.system("rm "+syspath+"unmapped"+str(k)+".txt") - else: - os.system("cut -f 1 "+syspath+"unmapped"+str(k)+".sam > "+syspath+"runmapped"+str(k)+".txt") - - os.system("rm "+syspath+"unmapped"+str(k)+".bam") - os.system("rm "+syspath+"umapped"+str(k)+".bam") - os.system("rm "+syspath+"rmapped"+str(k)+".bam") - os.system("python "+syspath+"seq_track.py "+syspath+"runmapped"+str(k)+".txt "+syspath+"seq"+str(k)+".fa "+syspath+"unmap_seq"+str(k)+".fa "+tp) #get unmapped sequence - os.system("python "+syspath+"truncate.py "+syspath+"unmap_seq"+str(k)+".fa "+shift+" "+syspath+"seq"+str(k+1)+".fa "+length) #truncate unmapped sequence - os.system("rm "+syspath+"seq"+str(k)+".fa") #Remove sequences being mapped - os.system("rm "+syspath+"map"+str(k)+".sam") #Remove mapping file - os.system("rm "+syspath+"unmap_seq"+str(k)+".fa") #Remove unmapped sequnce - os.system("rm "+syspath+"runmapped"+str(k)+".txt") - os.system("rm "+syspath+"unmapped"+str(k)+".sam") - - os.system("wc -l "+syspath+"seq"+str(k+1)+".fa > "+syspath+"count"+str(k+1)+".txt") - c = read_sp_file(syspath+"count"+str(k+1)+".txt") - if c[0][0] == '0': #If no reads is in the sequence file, stop - os.system("rm "+syspath+"count"+str(k+1)+".txt") - os.system("rm "+syspath+"seq"+str(k+1)+".fa") - break - os.system("rm "+syspath+"count"+str(k+1)+".txt") - k = k+1 - -ss = "" -for i in range(0,k+1): - ss = ss+" "+syspath+"mapped"+str(i)+".bam" - - -os.system("samtools merge -f "+output_file+" "+ss) -#print("samtools merge mapped_all.bam"+ss) -os.system("rm "+syspath+"mapped*.bam") -os.system("rm "+syspath+"ref*") - -