Mercurial > repos > bornea > saint_preprocessing
diff mzID_process2.py @ 72:5ec0b997fb13 draft
Uploaded
author | bornea |
---|---|
date | Sat, 27 Aug 2016 23:26:10 -0400 |
parents | f6d916d1d304 |
children | 50391fdc229a |
line wrap: on
line diff
--- a/mzID_process2.py Sat Aug 27 23:22:20 2016 -0400 +++ b/mzID_process2.py Sat Aug 27 23:26:10 2016 -0400 @@ -6,7 +6,7 @@ """ ####################################################################################### ## Description: ## -# This program will create inter, prey, and bait files from mzIdentML files +#This program will create inter, prey, and bait files from mzIdentML files ## Required input: ## # 1) mzIdentML file to be reformatted # 2) minimum PSM for quantification @@ -18,95 +18,95 @@ ins_path = sys.argv[5] class ReturnValue1(object): - def __init__(self, sequence, gene): - self.seqlength = sequence - self.genename = gene + def __init__(self, sequence, gene): + self.seqlength = sequence + self.genename = gene class ReturnValue2(object): - def __init__(self, inter, accessions): - self.inter = inter - self.accessions = accessions + def __init__(self, inter, accessions): + self.inter = inter + self.accessions = accessions def read_tab(infile): - with open(infile,'r') as x: - output = [] - for line in x: - line = line.strip() - temp = line.split('\t') - output.append(temp) - return output + with open(infile,'r') as x: + output = [] + for line in x: + line = line.strip() + temp = line.split('\t') + output.append(temp) + return output def printProgress (iteration, total, prefix = '', suffix = '', decimals = 1, barLength = 100): - """ - Call in a loop to create terminal progress bar - @params: - iteration - Required : current iteration (Int) - total - Required : total iterations (Int) - prefix - Optional : prefix string (Str) - suffix - Optional : suffix string (Str) - decimals - Optional : positive number of decimals in percent complete (Int) - barLength - Optional : character length of bar (Int) - """ - formatStr = "{0:." + str(decimals) + "f}" - percents = formatStr.format(100 * (iteration / float(total))) - filledLength = int(round(barLength * iteration / float(total))) - bar = '=' * filledLength + '-' * (barLength - filledLength) - sys.stdout.write('\r%s |%s| %s%s %s' % (prefix, bar, percents, '%', suffix)), - sys.stdout.flush() - if iteration == total: - sys.stdout.write('\n') - sys.stdout.flush() + """ + Call in a loop to create terminal progress bar + @params: + iteration - Required : current iteration (Int) + total - Required : total iterations (Int) + prefix - Optional : prefix string (Str) + suffix - Optional : suffix string (Str) + decimals - Optional : positive number of decimals in percent complete (Int) + barLength - Optional : character length of bar (Int) + """ + formatStr = "{0:." + str(decimals) + "f}" + percents = formatStr.format(100 * (iteration / float(total))) + filledLength = int(round(barLength * iteration / float(total))) + bar = '=' * filledLength + '-' * (barLength - filledLength) + sys.stdout.write('\r%s |%s| %s%s %s' % (prefix, bar, percents, '%', suffix)), + sys.stdout.flush() + if iteration == total: + sys.stdout.write('\n') + sys.stdout.flush() def get_info(uniprot_accession_in,fasta_db): - # Get aminoacid lengths and gene name. - error = open('error proteins.txt', 'a+') - data = open(fasta_db, 'r') - data_lines = data.readlines() - db_len = len(data_lines) - seqlength = 0 - count = 0 - last_line = data_lines[-1] - for data_line in data_lines: - if ">sp" in data_line: - namer = data_line.split("|")[2] - if uniprot_accession_in == data_line.split("|")[1]: - match = count + 1 - if 'GN=' in data_line: - lst = data_line.split('GN=') - lst2 = lst[1].split(' ') - genename = lst2[0] - if 'GN=' not in data_line: - genename = 'NA' - while ">sp" not in data_lines[match]: - if match <= db_len: - seqlength = seqlength + len(data_lines[match].strip()) - if data_lines[match] == last_line: - break - match = match + 1 - else: - break - return ReturnValue1(seqlength, genename) - if uniprot_accession_in == namer.split(" ")[0]: - match = count + 1 - # Ensures consistent spacing throughout. - if 'GN=' in data_line: - lst = data_line.split('GN=') - lst2 = lst[1].split(' ') - genename = lst2[0] - if 'GN=' not in data_line: - genename = 'NA' - while ">sp" not in data_lines[match]: - if match <= db_len: - seqlength = seqlength + len(data_lines[match].strip()) - if data_lines[match] == last_line: - break - match = match + 1 - else: - break - return ReturnValue1(seqlength, genename) - count = count + 1 - if seqlength == 0: - error.write(uniprot_accession_in + '\t' + "Uniprot not in Fasta" + '\n') - error.close - seqlength = 'NA' - genename = 'NA' - return ReturnValue1(seqlength, genename) + # Get aminoacid lengths and gene name. + error = open('error proteins.txt', 'a+') + data = open(fasta_db, 'r') + data_lines = data.readlines() + db_len = len(data_lines) + seqlength = 0 + count = 0 + last_line = data_lines[-1] + for data_line in data_lines: + if ">sp" in data_line: + namer = data_line.split("|")[2] + if uniprot_accession_in == data_line.split("|")[1]: + match = count + 1 + if 'GN=' in data_line: + lst = data_line.split('GN=') + lst2 = lst[1].split(' ') + genename = lst2[0] + if 'GN=' not in data_line: + genename = 'NA' + while ">sp" not in data_lines[match]: + if match <= db_len: + seqlength = seqlength + len(data_lines[match].strip()) + if data_lines[match] == last_line: + break + match = match + 1 + else: + break + return ReturnValue1(seqlength, genename) + if uniprot_accession_in == namer.split(" ")[0]: + match = count + 1 + # Ensures consistent spacing throughout. + if 'GN=' in data_line: + lst = data_line.split('GN=') + lst2 = lst[1].split(' ') + genename = lst2[0] + if 'GN=' not in data_line: + genename = 'NA' + while ">sp" not in data_lines[match]: + if match <= db_len: + seqlength = seqlength + len(data_lines[match].strip()) + if data_lines[match] == last_line: + break + match = match + 1 + else: + break + return ReturnValue1(seqlength, genename) + count = count + 1 + if seqlength == 0: + error.write(uniprot_accession_in + '\t' + "Uniprot not in Fasta" + '\n') + error.close + seqlength = 'NA' + genename = 'NA' + return ReturnValue1(seqlength, genename) def make_inter(mzIdentML,replicate,grouping): accession_index = mzIdentML[0].index("accession") PSMs = {} @@ -135,7 +135,7 @@ make_prey = sys.argv[3] db = sys.argv[4] if db == "None": - db = str(ins_path) + "/SwissProt_HUMAN_2015_12.fasta" + db = str(ins_path) + "/SwissProt_HUMAN_2015_12.fasta" make_bait = sys.argv[6] bait_bool = sys.argv[7] prey_file = sys.argv[8] @@ -143,37 +143,37 @@ inter_out = sys.argv[10] def bait_create(baits, infile): - # Verifies the Baits are valid in the Scaffold file and writes the Bait.txt. - baits = make_bait.split() - i = 0 - bait_file_tmp = open("bait.txt", "w") - order = [] - bait_cache = [] - while i < len(baits): - if baits[i+2] == "true": - T_C = "C" - else: - T_C = "T" - bait_line = baits[i] + "\t" + baits[i+1] + "\t" + T_C + "\n" - bait_cache.append(str(bait_line)) - i = i + 3 + # Verifies the Baits are valid in the Scaffold file and writes the Bait.txt. + baits = make_bait.split() + i = 0 + bait_file_tmp = open("bait.txt", "w") + order = [] + bait_cache = [] + while i < len(baits): + if baits[i+2] == "true": + T_C = "C" + else: + T_C = "T" + bait_line = baits[i] + "\t" + baits[i+1] + "\t" + T_C + "\n" + bait_cache.append(str(bait_line)) + i = i + 3 - for cache_line in bait_cache: - bait_file_tmp.write(cache_line) + for cache_line in bait_cache: + bait_file_tmp.write(cache_line) - bait_file_tmp.close() + bait_file_tmp.close() if bait_bool == 'false': - bait_create(make_bait, infile) - bait = "bait.txt" + bait_create(make_bait, infile) + bait = "bait.txt" else: - bait_temp_file = open(sys.argv[2], 'r') - bait_cache = bait_temp_file.readlines() - bait_file_tmp = open("bait.txt", "wr") - for cache_line in bait_cache: - bait_file_tmp.write(cache_line) - bait_file_tmp.close() - bait = "bait.txt" + bait_temp_file = open(sys.argv[2], 'r') + bait_cache = bait_temp_file.readlines() + bait_file_tmp = open("bait.txt", "wr") + for cache_line in bait_cache: + bait_file_tmp.write(cache_line) + bait_file_tmp.close() + bait = "bait.txt" bait = read_tab("bait.txt") inter = "" @@ -184,9 +184,9 @@ os.system(cmd) mzIdentML = read_tab("flat_mzIdentML.txt") inter = inter + make_inter(mzIdentML,bait[cnt][0],bait[cnt][1]).inter - print inter + print inter accessions.append(make_inter(mzIdentML,bait[cnt][0],bait[cnt][1]).accessions) - print accessions + print accessions cnt+=1 with open("inter.txt","w") as x: @@ -205,11 +205,11 @@ for i in unique_accessions: prey = prey + i + "\t" + str(get_info(i,db).seqlength) + "\t" + get_info(i,db).genename + "\n" start+=1 - printProgress(start, end) + printProgress(start, end) with open("prey.txt","w") as x: x.write(prey) os.rename("bait.txt", bait_out) os.rename("inter.txt", inter_out) if str(prey_file) != "None": - os.rename("prey.txt", prey_file) \ No newline at end of file + os.rename("prey.txt", prey_file) \ No newline at end of file