changeset 24:1d83dbbee373

Fixing a merge conflict
author melissacline
date Mon, 20 Jul 2015 15:41:22 -0700
parents 7346e5bab4ec (diff) d0674221a6ae (current diff)
children 80b6b64e8b5e
files mergeGenomicMatrixFiles.py
diffstat 4 files changed, 311 insertions(+), 24 deletions(-) [+]
line wrap: on
line diff
--- a/mergeGenomicFiles.xml	Wed Feb 11 17:28:52 2015 -0800
+++ b/mergeGenomicFiles.xml	Mon Jul 20 15:41:22 2015 -0700
@@ -1,20 +1,31 @@
-<tool id="mergeGenomicFiles" description="Merge two genomic datasets into a new dataset" name="mergeGenomicFiles" version="0.0.1">
+<tool id="mergeGenomicFiles" description="Merge two genomic datasets into a new dataset" name="Merge Genomic Datasets" version="0.0.1">
   <description>
-    Given two genomic datasets, merge them to create a third dataset with the row and column identifiers from both datasets.
+    Given two genomic datasets, merge them to create a larger dataset with the row and column identifiers from both datasets. Output this larger dataset, along with a 2-column matrix indicating the source file of each sample
   </description>
   <command interpreter="python">
-      mergeGenomicMatrixFiles.py $outputC $inputA $inputB
+      mergeGenomicMatrixFiles.py $inputA $inputB $outputC $outputSourceMatrix
+      #if $labelForDatasetA
+          --aLabel "${labelForDatasetA}"
+      #end if
+      #if $labelForDatasetB
+          --bLabel "${labelForDatasetB}"
+      #end if
   </command>
   <inputs>
     <param name="inputA" format="tabular" type="data" label="Genomic Dataset A"/>
-    <param name="inputB" format="tabular" type="data" label="Genomic Dataset B"/>
-  </inputs>
+    <param type="text" name="labelForDatasetA"  label="Dataset A Label (optional)" optional="true"/>
+    <param name="inputB" format="tabular" type="data" label="Genomic Dataset B"/> 
+    <param type="text" name="labelForDatasetB"  label="Dataset B Label (optional)" optional="true"/>
+ </inputs>
   <outputs>
-    <data name="outputC" format="tabular"/>
+    <data name="outputSourceMatrix" format="tabular" label="Genomic Data Sources"/>
+    <data name="outputC" format="tabular" label="Merged Genomic Data"/>
   </outputs>
   <help>
     ***Merge Genomic Datasets***
 
-    Given two genomic datasets, merge them to produce a third dataset that is the union of the first two.  The new dataset will contain all column labels from either dataset, and all row labels from either dataset.  If a row label appears in both datasets, the output dataset will contain, for that row, all values for the first set of columns, plus all values for the second set of columns.  If a row label appears in the first dataset only, the output dataset will contain the values for the columns of the first dataset, and blanks (indicating missing values) for the columns of the second da
+    Given two genomic datasets, merge them to produce a third dataset that is the union of the first two.  The new dataset will contain all column labels from either dataset, and all row labels from either dataset.  If a row label appears in both datasets, the output dataset will contain, for that row, all values for the first set of columns, plus all values for the second set of columns.  If a row label appears in the first dataset only, the output dataset will contain the values for the columns of the first dataset, and blanks (indicating missing values) for the columns of the second dataset.
+
+    To maintain provenance, this script also outputs a second matrix, with one row for each column in the output dataset, and two columns per row indicating which input dataset that column came from.  By default, the input dataset name is used to indicate which input file each column came from.  Optionally, the user can specify descriptive labels to be used in place of the filenames.  This all assumes that each column exists in only one input dataset.
   </help>
 </tool>
--- a/mergeGenomicMatrixFiles.py	Wed Feb 11 17:28:52 2015 -0800
+++ b/mergeGenomicMatrixFiles.py	Mon Jul 20 15:41:22 2015 -0700
@@ -1,21 +1,27 @@
 #!/usr/bin/env python
 
+import argparse
 import string,os,sys
 
-def header (samples, infile):
-    fin= open(infile,'r')
+def header (samples, sourceFiles, infile, labelThisFile):
+    if labelThisFile == None:
+        labelToUse = infile
+    else:
+        labelToUse = labelThisFile
+    fin= open(infile, 'U')
     #header, samples
     newSamples = string.split(string.strip(fin.readline()),'\t')[1:]
     for sample in newSamples:
         if sample not in samples:
             samples[sample]= len(samples)
+            sourceFiles[sample] = labelToUse
     fin.close()
     return
 
 def process(genes, samples, dataMatrix, infile):
     maxLength= len(samples)
 
-    fin= open(infile,'r')
+    fin= open(infile,'U')
     #header 
     newSamples = string.split(string.strip(fin.readline()),'\t')
     
@@ -41,7 +47,17 @@
     fin.close()
     return
 
-def outputMatrix(dataMatrix, samples, genes, outfile):
+
+def outputSourceMatrix(sourceData, outputFileName):
+    fout = open(outputFileName, "w")
+    fout.write("Sample\tSource\n")
+    for thisSample in sourceData.keys():
+        fout.write("%s\t%s\n" % (thisSample, sourceData[thisSample]))
+    fout.close()
+    return
+
+
+def outputMergedMatrix(dataMatrix, samples, genes, outfile):
     fout = open(outfile,"w")
     maxLength= len(samples)
     sList=[]
@@ -66,23 +82,38 @@
     return
 
 if __name__ == '__main__' :
-    if len(sys.argv[:]) <4:
-        print "python mergeFilesByColumn.py output inputfile(s)"
-        print "**********memory intensive, not for very genomic data with hugo number of probes"
-        print "this is merging data A+B=C\n"
-        sys.exit()
-
-    inFiles = sys.argv[2:]
-    outfile = sys.argv[1]
+    #
+    # The input files to this script are two or more matrices, in which
+    # columns represent samples and rows represent genes or measurements.
+    # There are two output files: outMergedData contains the input data merged
+    # into a single matrix, and outSourceMatrix is a two-column matrix 
+    # indicating which file each sample (or column label) came from.  This
+    # assumes that each sample came from at most one file.
+    #
+    parser = argparse.ArgumentParser()
+    parser.add_argument("inFileA", type=str, help="First input file")
+    parser.add_argument("inFileB", type=str, help="Second input file")
+    parser.add_argument("outMergedData", type=str, 
+                        help="Filename for the merged dataset")
+    parser.add_argument("outSourceMatrix", type=str,
+                        help="""Filename for a Nx2 matrix that indicates
+                                the source file of each column""")
+    parser.add_argument("--aLabel", type=str, default=None,
+                        help="User-friendly label for the first input file")
+    parser.add_argument("--bLabel", type=str, default=None,
+                        help="User-friendly label for the second input file")
+    args = parser.parse_args()
 
     genes={}
     samples={}
+    sourceFiles = {}
     dataMatrix=[]
 
-    for infile in inFiles:
-        header (samples, infile)
+    header(samples, sourceFiles, args.inFileA, args.aLabel)
+    header(samples, sourceFiles, args.inFileB, args.bLabel)
 
-    for infile in inFiles:
-        process(genes, samples, dataMatrix, infile)
+    process(genes, samples, dataMatrix, args.inFileA)
+    process(genes, samples, dataMatrix, args.inFileB)
 
-    outputMatrix(dataMatrix, samples, genes, outfile)
+    outputSourceMatrix(sourceFiles, args.outSourceMatrix)
+    outputMergedMatrix(dataMatrix, samples, genes, args.outMergedData)
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mergeMutationDatasets.xml	Mon Jul 20 15:41:22 2015 -0700
@@ -0,0 +1,31 @@
+<tool id="mergeMutationDatasets" description="Merge two Xena positional mutation datasets into a new dataset" name="Merge Positional Mutation Data" version="0.0.1">
+  <description>
+    Given two mutation datasets, merge them to create a larger dataset with the mutations from both datasets. Output this larger dataset, along with a 2-column matrix indicating the source of each mutation
+  </description>
+  <command interpreter="python">
+      mergeXenaMutation.py $outputC $outputSourceMatrix $errorLog  $inputA $inputB 
+      #if $labelForDatasetA
+          --aLabel "${labelForDatasetA}"
+      #end if
+      #if $labelForDatasetB
+          --bLabel "${labelForDatasetB}"
+      #end if
+  </command>
+  <inputs>
+    <param name="inputA" format="tabular" type="data" label="Mutation Dataset A"/>
+    <param type="text" name="labelForDatasetA"  label="Dataset A Label (optional)" optional="true"/>
+    <param name="inputB" format="tabular" type="data" label="Mutation Dataset B"/> 
+    <param type="text" name="labelForDatasetB"  label="Dataset B Label (optional)" optional="true"/>
+ </inputs>
+  <outputs>
+    <data name="errorLog" format="data" label="Execution Log"/>
+    <data name="outputSourceMatrix" format="tabular" label="Mutation Data Sources"/> 
+    <data name="outputC" format="tabular" label="Merged Mutation Data"/>
+  </outputs>
+  <help>
+    ***Merge Xena Positional Mutation Datasets***
+
+    Given two datasets of mutation data as formatted for the UCSC Xena Browser, merge them to produce a third dataset that is the union of the first two.  The new dataset will contain all mutations from either dataset. 
+
+    To maintain provenance, this script also outputs a second matrix, with one row for each sample ID that appears in the output dataset, and two columns per row indicating which input dataset(s) contained some mutation data for that sample.  By default, the input dataset name is used to indicate which input file each column came from.  Optionally, the user can specify descriptive labels to be used in place of the dataset names.   </help>
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mergeXenaMutation.py	Mon Jul 20 15:41:22 2015 -0700
@@ -0,0 +1,214 @@
+#!/usr/bin/env python
+
+import argparse
+import string, os, sys
+
+requiredCOLs = ["chr", "start","end","reference","alt","gene","effect"]
+
+def headerError(filename, column, ferror):
+  ferror.write(filename +" does not have column " + column+"\n")
+  ferror.close()
+  sys.exit(1)
+
+def findAnyValueInList (values, dataList):
+  for value in values:
+    for i in range(0,len(dataList)):
+      if value == dataList[i]:
+        return i
+  return -1
+
+def header (infile, ferror):
+  fin= open(infile,'U')
+
+  columnDic ={}
+  #header
+  line = fin.readline()
+  fin.close()
+  if line [0]=="#":
+    line = line[1:-1]
+  data = string.split(line,"\t")
+
+  columnDic["chr"]= findAnyValueInList (["chr","chrom"], data)
+  if columnDic["chr"] ==-1:
+    headerError(infile, "chr", ferror)
+
+  columnDic["start"]= findAnyValueInList (["start","chrStart"], data)
+  if columnDic["start"] == -1:
+    headerError(infile, "start", ferror)
+
+  columnDic["end"]= findAnyValueInList (["end","chrEnd"], data)
+  if columnDic["end"] == -1:
+    headerError(infile, "end", ferror)
+
+  columnDic["alt"]= findAnyValueInList (["alt"], data)
+  if columnDic["alt"] == -1:
+    headerError(infile, "alt", ferror)
+
+  columnDic["reference"]= findAnyValueInList (["reference","ref"], data)
+  if columnDic["reference"] == -1:
+    headerError(infile, "reference", ferror)
+
+  columnDic["gene"]= findAnyValueInList (["gene"], data)
+  if columnDic["gene"] == -1:
+    headerError(infile, "gene", ferror)
+
+  columnDic["effect"]= findAnyValueInList (["effect"], data)
+  if columnDic["effect"] == -1:
+    headerError(infile, "effect", ferror)
+
+  requiredCols = columnDic.keys()
+  requiredColsPos = columnDic.values()
+  for i in range(1,len(data)):
+    if i not in requiredColsPos:
+      columnDic [data[i]]=i
+
+  return columnDic
+
+def summarizeColumns(inFiles, fileColumn, allCols, ferror):
+  for infile in inFiles:
+    columnDic = header (infile, ferror)
+    fileColumn [infile] = columnDic
+    for col in columnDic:
+      if col not in allCols:
+        allCols.append(col)
+  return
+
+def outputHeader (requiredCOLs,allCols,fout):
+  fout.write("#sample")
+  for col in requiredCOLs:
+    fout.write("\t"+col)
+  for col in allCols:
+    if col not in requiredCOLs:
+      fout.write("\t"+col)
+  fout.write("\n")
+  fout.close()
+  return
+
+def processAndOutput(infile,requiredCOLs,allCols,columnDic,fout):
+  fin = open(infile,'U')
+  fin.readline()
+  while 1:
+    line = fin.readline()[:-1]
+    if line =="":
+      break
+    data = string.split(line,'\t')
+    fout.write(data[0])
+    for col in requiredCOLs:
+      pos = columnDic[col]
+      fout.write("\t"+ data[pos])
+    for col in allCols:
+      if col not in requiredCOLs:
+        if col in columnDic:
+          pos = columnDic[col]
+          fout.write("\t"+ data[pos])
+        else:
+          fout.write("\t")
+    fout.write("\n")
+  fin.close()
+  return
+
+def collectSource(inFile, label, sampleDic):
+  fin = open(inFile,'U')
+  fin.readline()
+  while 1:
+    line = fin.readline()[:-1]
+    if line =="":
+      break
+    sample = string.split(line,'\t')[0]
+    if sample not in sampleDic:
+      sampleDic[sample]=[]
+    if inFile not in sampleDic[sample]:
+      sampleDic[sample].append(label)
+  fin.close()
+  return
+
+def outputSampleDic (sampleDic, outPhenotypeFile):
+  fout = open(outPhenotypeFile,'w')
+  fout.write("sample\tsource\n")
+  for sample in sampleDic:
+    source = sampleDic[sample]
+    source.sort()
+    fout.write(sample+"\t"+string.join(source,", ")+"\n")
+  fout.close()
+  return
+
+if __name__ == '__main__' :
+  if len(sys.argv[:]) <6:
+    print "python mergeMultipleXenaMutation.py outputXenaMutation outputPhenotypeMatrix errorLog inputfile(s)"
+    print "this is merging data A+B=C for mutation by position type of data\n"
+    sys.exit(1)
+
+  #
+  # The input files to this script are two or more matrices, in which
+  # columns represent samples and rows represent genes or measurements.
+  # There are two output files: outMergedData contains the input data merged
+  # into a single matrix, and outSourceMatrix is a two-column matrix
+  # indicating which file each sample (or column label) came from.  This
+  # assumes that each sample came from at most one file.
+  #
+  parser = argparse.ArgumentParser()
+  parser.add_argument("outMergedData", type=str,
+                      help="Filename for the merged dataset")
+  parser.add_argument("outSourceMatrix", type=str,
+                      help="""Filename for a Nx2 matrix that indicates
+                                the source file of each column""")
+  parser.add_argument("errorLog", type=str,
+                      help="""Error log""")
+  parser.add_argument("inFileA", type=str, help="First input file")
+  parser.add_argument("inFileB", type=str, help="Second input file")
+  parser.add_argument("--aLabel", type=str, default=None,
+                      help="User-friendly label for the first input file")
+  parser.add_argument("--bLabel", type=str, default=None,
+                      help="User-friendly label for the second input file")
+  args = parser.parse_args()
+
+                                    
+  #inFiles = sys.argv[4:]
+  inFiles = list()
+  inFiles.append(args.inFileA)
+  inFiles.append(args.inFileB)
+  errofile = args.errorLog
+  outfile = args.outMergedData
+  #print outfile
+  outPhenotypeFile = args.outSourceMatrix
+  #print outPhenotypeFile
+    
+  ferror = open(errofile,'w')
+    
+  #get all the columns, build fileColumn dictionary
+  fileColumn={}
+  allCols =[]
+  summarizeColumns(inFiles, fileColumn, allCols, ferror)
+  ferror.close()
+
+  #output header line
+  fout = open(outfile,'w')
+  outputHeader (requiredCOLs,allCols,fout)
+  
+  #process and output combined mutationXena file
+  fout = open(outfile,'a')
+
+  columnDic = fileColumn[args.inFileA]
+  processAndOutput(args.inFileA,requiredCOLs,allCols,columnDic,fout)
+  columnDic = fileColumn[args.inFileB]
+  processAndOutput(args.inFileB,requiredCOLs,allCols,columnDic,fout)
+  fout.close()
+  
+  #collect sample from source information
+  sampleDic ={}
+  if args.aLabel is None:
+    collectSource(args.inFileA, args.inFileA, sampleDic)
+  else:
+    collectSource(args.inFileA, args.aLabel, sampleDic)
+  if args.bLabel is None:
+    collectSource(args.inFileB, args.inFileB, sampleDic)
+  else:
+    collectSource(args.inFileB, args.bLabel, sampleDic)
+
+
+  #output sample source information as phenotype matrix
+  outputSampleDic (sampleDic, outPhenotypeFile)
+
+
+
+