Mercurial > repos > sagun98 > micropita
diff galaxy_micropita/src/breadcrumbs/scripts/scriptPcoa.py @ 3:8fb4630ab314 draft default tip
Uploaded
author | sagun98 |
---|---|
date | Thu, 03 Jun 2021 17:07:36 +0000 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/galaxy_micropita/src/breadcrumbs/scripts/scriptPcoa.py Thu Jun 03 17:07:36 2021 +0000 @@ -0,0 +1,143 @@ +#!/usr/bin/env python +""" +Author: Timothy Tickle +Description: Make PCoA of an abundance file +""" + +__author__ = "Timothy Tickle" +__copyright__ = "Copyright 2012" +__credits__ = ["Timothy Tickle"] +__license__ = "" +__version__ = "" +__maintainer__ = "Timothy Tickle" +__email__ = "ttickle@sph.harvard.edu" +__status__ = "Development" + +import sys +import argparse +from src.AbundanceTable import AbundanceTable +from src.Metric import Metric +import csv +import os +from src.PCoA import PCoA + +#Set up arguments reader +argp = argparse.ArgumentParser( prog = "scriptPcoa.py", + description = """PCoAs an abundance file given a metadata.\nExample:python scriptPcoa.py -i TID -l STSite""" ) + +#Arguments +#For table +argp.add_argument("-i","--id", dest="sIDName", default="ID", help="Abundance Table ID") +argp.add_argument("-l","--meta", dest="sLastMetadataName", help="Last metadata name") +argp.add_argument("-d","--fDelim", dest= "cFileDelimiter", action= "store", default="\t", help="File delimiter, default tab") +argp.add_argument("-f","--featureDelim", dest="cFeatureNameDelimiter", action= "store", metavar="Feature Name Delimiter", default="|", help="Feature delimiter") + +argp.add_argument("-n","--doNorm", dest="fDoNormData", action="store_true", default=False, help="Flag to turn on normalization") +argp.add_argument("-s","--doSum", dest="fDoSumData", action="store_true", default=False, help="Flag to turn on summation") + +argp.add_argument("-p","--paint", dest="sLabel", metavar= "Label", default=None, help="Label to paint in the PCoA") +argp.add_argument("-m","--metric", dest="strMetric", metavar = "distance", default = PCoA.c_BRAY_CURTIS, help ="Distance metric to use. Pick from braycurtis, canberra, chebyshev, cityblock, correlation, cosine, euclidean, hamming, spearman, sqeuclidean, unifrac_unweighted, unifrac_weighted") +argp.add_argument("-o","--outputFile", dest="strOutFile", metavar= "outputFile", default=None, help="Specify the path for the output figure.") +argp.add_argument("-D","--DistanceMatrix", dest="strFileDistanceMatrix", metavar= "strFileDistanceMatrix", default=None, help="Specify the path for outputing the distance matrix (if interested). Default this will not output.") +argp.add_argument("-C","--CoordinatesMatrix", dest="strFileCoordinatesMatrix", metavar= "strFileCoordinatesMatrix", default=None, help="Specify the path for outputing the x,y coordinates matrix (Dim 1 and 2). Default this will not output.") + +# Unifrac arguments +argp.add_argument("-t","--unifracTree", dest="istrmTree", metavar="UnifracTreeFile", default=None, help="Optional file only needed for UniFrac calculations.") +argp.add_argument("-e","--unifracEnv", dest="istrmEnvr", metavar="UnifracEnvFile", default=None, help="Optional file only needed for UniFrac calculations.") +argp.add_argument("-c","--unifracColor", dest="fileUnifracColor", metavar="UnifracColorFile", default = None, help="A text file indicating the groupings of metadata to color. Each line in the file is a group to color. An example file line would be 'GroupName:ID,ID,ID,ID'") + +argp.add_argument("strFileAbund", metavar = "Abundance file", nargs="?", help ="Input data file") + +args = argp.parse_args( ) + +#Read in abundance table +abndTable = None +if args.strFileAbund: + abndTable = AbundanceTable.funcMakeFromFile(args.strFileAbund, + cDelimiter = args.cFileDelimiter, + sMetadataID = args.sIDName, + sLastMetadata = args.sLastMetadataName, + cFeatureNameDelimiter= args.cFeatureNameDelimiter) + + #Normalize if need + if args.fDoSumData: + abndTable.funcSumClades() + + #Sum if needed + if args.fDoNormData: + abndTable.funcNormalize() + +#Get the metadata to paint +lsKeys = None +if abndTable: + lsKeys = abndTable.funcGetMetadataCopy().keys() if not args.sLabel else [args.sLabel] + +#Get pieces of output file +if not args.strOutFile: + if not args.strFileAbund: + args.strOutFile = os.path.splitext(os.path.basename(args.istrmEnvr))[0]+"-pcoa.pdf" + else: + args.strOutFile = os.path.splitext(os.path.basename(args.strFileAbund))[0]+"-pcoa.pdf" +lsFilePieces = os.path.splitext(args.strOutFile) + +# Make PCoA object +# Get PCoA object and plot +pcoa = PCoA() +if(not args.strMetric in [Metric.c_strUnifracUnweighted,Metric.c_strUnifracWeighted]) and abndTable: + pcoa.loadData(abndTable,True) +# Optional args.strFileDistanceMatrix if not none will force a printing of the distance measures to the path in args.strFileDistanceMatrix +pcoa.run(tempDistanceMetric=args.strMetric, iDims=2, strDistanceMatrixFile=args.strFileDistanceMatrix, istrmTree=args.istrmTree, istrmEnvr=args.istrmEnvr) + +# Write dim 1 and 2 coordinates to file +if args.strFileCoordinatesMatrix: + lsIds = pcoa.funcGetIDs() + mtrxCoordinates = pcoa.funcGetCoordinates() + csvrCoordinates = csv.writer(open(args.strFileCoordinatesMatrix, 'w')) + csvrCoordinates.writerow(["ID","Dimension_1","Dimension_2"]) + for x in xrange(mtrxCoordinates.shape[0]): + strId = lsIds[x] if lsIds else "" + csvrCoordinates.writerow([strId]+mtrxCoordinates[x].tolist()) + +# Paint metadata +if lsKeys: + for iIndex in xrange(len(lsKeys)): + lsMetadata = abndTable.funcGetMetadata(lsKeys[iIndex]) + + pcoa.plotList(lsLabelList = lsMetadata, + strOutputFileName = lsFilePieces[0]+"-"+lsKeys[iIndex]+lsFilePieces[1], + iSize=20, + dAlpha=1.0, + charForceColor=None, + charForceShape=None, + fInvert=False, + iDim1=1, + iDim2=2) + +if args.strMetric in [Metric.c_strUnifracUnweighted,Metric.c_strUnifracWeighted]: + + c_sNotGiven = "Not_specified" + + lsIds = pcoa.funcGetIDs() + lsGroupLabels = [c_sNotGiven for s in lsIds] + + if args.fileUnifracColor: + + # Read color file and make a dictionary to convert ids + lsColorLines = csv.reader(open(args.fileUnifracColor)) + dictConvertIDToGroup = {} + for lsLine in lsColorLines: + if lsLine: + sGroupID, sFirstID = lsLine[0].split(":") + dictConvertIDToGroup.update(dict([(sID,sGroupID) for sID in [sFirstID]+lsLine[1:]])) + + lsGroupLabels = [dictConvertIDToGroup.get(sID,c_sNotGiven) for sID in lsIds] + + pcoa.plotList(lsLabelList = lsGroupLabels, + strOutputFileName = lsFilePieces[0]+"-"+args.strMetric+lsFilePieces[1], + iSize=20, + dAlpha=1.0, + charForceColor=None, + charForceShape=None, + fInvert=False, + iDim1=1, + iDim2=2)