Mercurial > repos > davidvanzessen > argalaxy_tools
view imgt_loader/imgt_loader.py @ 1:90a05ff900db draft
Uploaded
author | davidvanzessen |
---|---|
date | Mon, 05 Sep 2016 04:46:06 -0400 |
parents | afe85eb6572e |
children |
line wrap: on
line source
import pandas as pd try: pd.options.mode.chained_assignment = None # default='warn' except: pass import re import argparse import os def stop_err( msg, ret=1 ): sys.stderr.write( msg ) sys.exit( ret ) #docs.python.org/dev/library/argparse.html parser = argparse.ArgumentParser() parser.add_argument("--summ", help="The 1_Summary file from the imgt output") parser.add_argument("--aa", help="The 5_AA-Sequence file from the imgt output") parser.add_argument("--junction", help="The 6_Junction file from the imgt output") parser.add_argument("--output", help="Output file") args = parser.parse_args() old_summary_columns = [u'Sequence ID', u'JUNCTION frame', u'V-GENE and allele', u'D-GENE and allele', u'J-GENE and allele', u'CDR1-IMGT length', u'CDR2-IMGT length', u'CDR3-IMGT length', u'Orientation'] old_sequence_columns = [u'CDR1-IMGT', u'CDR2-IMGT', u'CDR3-IMGT'] old_junction_columns = [u'JUNCTION'] added_summary_columns = [u'Functionality', u'V-REGION identity %', u'V-REGION identity nt', u'D-REGION reading frame', u'AA JUNCTION', u'Functionality comment', u'Sequence'] added_sequence_columns = [u'FR1-IMGT', u'FR2-IMGT', u'FR3-IMGT', u'CDR3-IMGT', u'JUNCTION', u'J-REGION', u'FR4-IMGT'] added_junction_columns = [u"P3'V-nt nb", u'N-REGION-nt nb', u'N1-REGION-nt nb', u"P5'D-nt nb", u"P3'D-nt nb", u'N2-REGION-nt nb', u"P5'J-nt nb", u"3'V-REGION trimmed-nt nb", u"5'D-REGION trimmed-nt nb", u"3'D-REGION trimmed-nt nb", u"5'J-REGION trimmed-nt nb", u"N-REGION", u"N1-REGION", u"N2-REGION"] outFile = args.output #fSummary = pd.read_csv(triplets[0][0], sep="\t", low_memory=False) fSummary = pd.read_csv(args.summ, sep="\t", dtype=object) #fSequence = pd.read_csv(triplets[0][1], sep="\t", low_memory=False) fSequence = pd.read_csv(args.aa, sep="\t", dtype=object) #fJunction = pd.read_csv(triplets[0][2], sep="\t", low_memory=False) fJunction = pd.read_csv(args.junction, sep="\t", dtype=object) tmp = fSummary[["Sequence ID", "JUNCTION frame", "V-GENE and allele", "D-GENE and allele", "J-GENE and allele"]] tmp["CDR1 Seq"] = fSequence["CDR1-IMGT"] tmp["CDR1 Length"] = fSummary["CDR1-IMGT length"] tmp["CDR2 Seq"] = fSequence["CDR2-IMGT"] tmp["CDR2 Length"] = fSummary["CDR2-IMGT length"] tmp["CDR3 Seq"] = fSequence["CDR3-IMGT"] tmp["CDR3 Length"] = fSummary["CDR3-IMGT length"] tmp["CDR3 Seq DNA"] = fJunction["JUNCTION"] tmp["CDR3 Length DNA"] = '1' tmp["Strand"] = fSummary["Orientation"] tmp["CDR3 Found How"] = 'a' for col in added_summary_columns: tmp[col] = fSummary[col] for col in added_sequence_columns: tmp[col] = fSequence[col] for col in added_junction_columns: tmp[col] = fJunction[col] outFrame = tmp outFrame.columns = [u'ID', u'VDJ Frame', u'Top V Gene', u'Top D Gene', u'Top J Gene', u'CDR1 Seq', u'CDR1 Length', u'CDR2 Seq', u'CDR2 Length', u'CDR3 Seq', u'CDR3 Length', u'CDR3 Seq DNA', u'CDR3 Length DNA', u'Strand', u'CDR3 Found How', u'Functionality', 'V-REGION identity %', 'V-REGION identity nt', 'D-REGION reading frame', 'AA JUNCTION', 'Functionality comment', 'Sequence', 'FR1-IMGT', 'FR2-IMGT', 'FR3-IMGT', 'CDR3-IMGT', 'JUNCTION', 'J-REGION', 'FR4-IMGT', 'P3V-nt nb', 'N-REGION-nt nb', 'N1-REGION-nt nb', 'P5D-nt nb', 'P3D-nt nb', 'N2-REGION-nt nb', 'P5J-nt nb', '3V-REGION trimmed-nt nb', '5D-REGION trimmed-nt nb', '3D-REGION trimmed-nt nb', '5J-REGION trimmed-nt nb', "N-REGION", "N1-REGION", "N2-REGION"] """ IGHV[0-9]-[0-9ab]+-?[0-9]?D? TRBV[0-9]{1,2}-?[0-9]?-?[123]? IGKV[0-3]D?-[0-9]{1,2} IGLV[0-9]-[0-9]{1,2} TRAV[0-9]{1,2}(-[1-46])?(/DV[45678])? TRGV[234589] TRDV[1-3] IGHD[0-9]-[0-9ab]+ TRBD[12] TRDD[1-3] IGHJ[1-6] TRBJ[12]-[1-7] IGKJ[1-5] IGLJ[12367] TRAJ[0-9]{1,2} TRGJP?[12] TRDJ[1-4] """ vPattern = [r"(IGHV[0-9]-[0-9ab]+-?[0-9]?D?)", r"(TRBV[0-9]{1,2}-?[0-9]?-?[123]?)", r"(IGKV[0-3]D?-[0-9]{1,2})", r"(IGLV[0-9]-[0-9]{1,2})", r"(TRAV[0-9]{1,2}(-[1-46])?(/DV[45678])?)", r"(TRGV[234589])", r"(TRDV[1-3])", r"(IGHV[0-9]S[0-9]+)"] dPattern = [r"(IGHD[0-9]-[0-9ab]+)", r"(TRBD[12])", r"(TRDD[1-3])"] jPattern = [r"(IGHJ[1-6])", r"(TRBJ[12]-[1-7])", r"(IGKJ[1-5])", r"(IGLJ[12367])", r"(TRAJ[0-9]{1,2})", r"(TRGJP?[12])", r"(TRDJ[1-4])"] vPattern = re.compile(r"|".join(vPattern)) dPattern = re.compile(r"|".join(dPattern)) jPattern = re.compile(r"|".join(jPattern)) def filterGenes(s, pattern): if type(s) is not str: return "NA" res = pattern.search(s) if res: return res.group(0) return "NA" outFrame["Top V Gene"] = outFrame["Top V Gene"].apply(lambda x: filterGenes(x, vPattern)) outFrame["Top D Gene"] = outFrame["Top D Gene"].apply(lambda x: filterGenes(x, dPattern)) outFrame["Top J Gene"] = outFrame["Top J Gene"].apply(lambda x: filterGenes(x, jPattern)) tmp = outFrame["VDJ Frame"] tmp = tmp.replace("in-frame", "In-frame") tmp = tmp.replace("null", "Out-of-frame") tmp = tmp.replace("out-of-frame", "Out-of-frame") outFrame["VDJ Frame"] = tmp outFrame["CDR3 Length DNA"] = outFrame["CDR3 Seq DNA"].map(str).map(len) safeLength = lambda x: len(x) if type(x) == str else 0 #outFrame = outFrame[(outFrame["CDR3 Seq DNA"].map(safeLength) > 0) & (outFrame["Top V Gene"] != "NA") & (outFrame["Top J Gene"] != "NA")] #filter out weird rows? #outFrame = outFrame[(outFrame["CDR3 Seq DNA"].map(safeLength) > 0) & (outFrame["Top V Gene"] != "NA") & (outFrame["Top D Gene"] != "NA") & (outFrame["Top J Gene"] != "NA")] #filter out weird rows? outFrame.to_csv(outFile, sep="\t", index=False, index_label="index")