annotate src/PCLToGraphlanCoreGene.py @ 8:e9677425c6c3 default tip

Updated the structure of the libraries
author george.weingart@gmail.com
date Mon, 09 Feb 2015 12:17:40 -0500
parents e0b5980139d9
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
1 #!/usr/bin/env python
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
2 #####################################################################################
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
3 #Copyright (C) <2012>
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
4 #
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
5 #Permission is hereby granted, free of charge, to any person obtaining a copy of
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
6 #this software and associated documentation files (the "Software"), to deal in the
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
7 #Software without restriction, including without limitation the rights to use, copy,
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
8 #modify, merge, publish, distribute, sublicense, and/or sell copies of the Software,
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
9 #and to permit persons to whom the Software is furnished to do so, subject to
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
10 #the following conditions:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
11 #
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
12 #The above copyright notice and this permission notice shall be included in all copies
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
13 #or substantial portions of the Software.
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
14 #
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
15 #THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED,
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
16 #INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
17 #PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
18 #HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
19 #OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
20 #SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
21 #
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
22 # This file is a component of the MaAsLin (Multivariate Associations Using Linear Models),
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
23 # authored by the Huttenhower lab at the Harvard School of Public Health
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
24 # (contact Timothy Tickle, ttickle@hsph.harvard.edu).
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
25 #####################################################################################
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
26
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
27 __author__ = "Timothy Tickle"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
28 __copyright__ = "Copyright 2012"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
29 __credits__ = ["Timothy Tickle"]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
30 __license__ = ""
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
31 __version__ = ""
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
32 __maintainer__ = "Timothy Tickle"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
33 __email__ = "ttickle@sph.harvard.edu"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
34 __status__ = "Development"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
35
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
36 import argparse
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
37 import csv
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
38 from operator import itemgetter
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
39 import re
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
40 import sys
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
41
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
42 #Helper function which returns a boolean indicator of an input string being parsable as an int
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
43 def funcIsInt(strInt):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
44 try:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
45 int(strInt)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
46 return True
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
47 except:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
48 return False
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
49
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
50 #Helper function that gets the index of the name and gives the last value of the list for - or the first value depending on the position
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
51 # This supports the ranging in the read.config files
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
52 #If no range is given then the result is just one index of the given name
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
53 def funcGetIndices(lsFeature, lsFunctionNames):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
54 if(len(lsFeature)) == 1:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
55 if(funcIsInt(lsFeature[0])):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
56 return int(lsFeature[0])-1
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
57 return [lsFeatureNames.index(lsFeature[0])]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
58 if(len(lsFeature)) == 2:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
59 iIndices = []
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
60 iPosition = 1
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
61 for sFeature in lsFeature:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
62 if(sFeature==""):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
63 if(iPosition==1):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
64 iIndices.append(2)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
65 elif(iPosition==2):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
66 iIndices.append(len(lsFunctionNames)-1)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
67 elif(funcIsInt(sFeature)):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
68 iIndices.append(int(sFeature)-1)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
69 else:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
70 iIndices.append(lsFeatureNames.index(sFeature))
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
71 iPosition = iPosition + 1
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
72 return iIndices
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
73
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
74 #Constants
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
75 #The line indicating the rows to read
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
76 c_MatrixName = "Matrix:"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
77 c_DataMatrix = "Abundance"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
78 c_strRows = "Read_PCL_Rows:"
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
79
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
80 #Set up arguments reader
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
81 argp = argparse.ArgumentParser( prog = "PCLToGraphlanCoreGene.py",
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
82 description = """Converts PCL files to Graphlan core gene files.""" )
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
83
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
84 #Arguments
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
85 argp.add_argument("strInputPCL", metavar = "PCLFile", type = argparse.FileType("r"), help ="Input PCl file used in maaslin")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
86 argp.add_argument("strInputRC", metavar = "RCFile", type = argparse.FileType("r"), help ="Input read config file used in maaslin")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
87 argp.add_argument("strOutputCoreGene", metavar = "CoreGeneFile", type = argparse.FileType("w"), help ="Output core gene file for graphlan")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
88
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
89 args = argp.parse_args( )
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
90
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
91 #Read in read config table and get the rows/columns to use
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
92 #Indicates if we are reading a data matrix
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
93 fIsData = False
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
94 #Holds the indices ranges
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
95 #List of lists,each internal list hold 1 or 2 indices, if two it indicates a range from the first to the second
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
96 llsIndices = []
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
97 csvRC = open(args.strInputRC,'r') if isinstance(args.strInputRC, str) else args.strInputRC
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
98 fRC = csv.reader(csvRC, delimiter=" ")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
99 for sLine in fRC:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
100 #Get the row indices or names
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
101 if len(sLine):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
102 if sLine[0] == c_MatrixName:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
103 fIsData = sLine[1] == c_DataMatrix
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
104 if sLine[0] == c_strRows:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
105 if fIsData:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
106 llsIndices = [sIndexRange.split("-") for sIndexRange in sLine[1].split(",")]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
107 break
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
108 csvRC.close()
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
109
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
110 # Check to make sure RC file is read
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
111 if len(llsIndices)==0:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
112 print("PCLToGraphlanCoreGene:: Could Not find indices in RC file "+args.strInputRC+".")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
113
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
114 #Read in the PCL file and parse the file names to core genes format
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
115 csvPCL = open(args.strInputPCL,'r') if isinstance(args.strInputPCL, str) else args.strInputPCL
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
116 fPCL = csv.reader(csvPCL,delimiter="\t")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
117 #The first column of the csv file
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
118 lsFeatureNames = [sLine[0] for sLine in fPCL]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
119 csvPCL.close()
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
120
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
121 # Check to make sure PCL file is read
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
122 if len(lsFeatureNames)==0:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
123 print("PCLToGraphlanCoreGene:: Could Not find features in PCL file "+args.strInputPCL+".")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
124
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
125 #If the indices are names switch with numbers otherwise subtract 1 because they are ment for R
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
126 liConvertedRangedIndices = [funcGetIndices(sIndex,lsFeatureNames) for sIndex in llsIndices] if len(llsIndices)>0 else []
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
127 llsIndices = None
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
128
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
129 #If there are any ranges, reduce to lists of indices
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
130 liConvertedIndices = []
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
131 for lsIndices in liConvertedRangedIndices:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
132 lsIndices.sort()
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
133 iLenIndices = len(lsIndices)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
134 if iLenIndices > 2:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
135 print "Error, received more than 2 indices in a range. Stopped."
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
136 exit()
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
137 liConvertedIndices.extend(lsIndices if iLenIndices == 1 else range(lsIndices[0],lsIndices[1]+1))
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
138 liConvertedRangedIndices = None
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
139
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
140 #Collapse all indices to a set which is then sorted
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
141 liConvertedIndices = sorted(list(set(liConvertedIndices)))
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
142
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
143 #Reduce name of features to just bugs indicated by indices
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
144 lsFeatureNames = itemgetter(*liConvertedIndices)(lsFeatureNames)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
145 liConvertedIndices = None
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
146
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
147 #Change the bug names to the correct formatting (clades seperated by .)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
148 lsFeatureNames = sorted(lsFeatureNames)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
149 lsFeatureNames = [re.sub("^[A-Za-z]__","",sBug) for sBug in lsFeatureNames]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
150 lsFeatureNames = [[re.sub("\|*[A-Za-z]__|\|",".",sBug)] for sBug in lsFeatureNames]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
151
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
152 #If this is an OTU, append the number and the genus level together for a more descriptive termal name
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
153 lsFeatureNamesModForOTU = []
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
154 for sBug in lsFeatureNames:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
155 lsBug = sBug[0].split(".")
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
156 if(len(lsBug))> 1:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
157 if(lsBug[-1].isdigit()):
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
158 lsBug[-2]=lsBug[-2]+"_"+lsBug[-1]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
159 lsBug = lsBug[0:-1]
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
160 lsFeatureNamesModForOTU.append([".".join(lsBug)])
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
161 else:
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
162 lsFeatureNamesModForOTU.append([lsBug[0]])
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
163
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
164 #Output core gene file
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
165 csvCG = open(args.strOutputCoreGene,'w') if isinstance(args.strOutputCoreGene, str) else args.strOutputCoreGene
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
166 fCG = csv.writer(csvCG)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
167 fCG.writerows(lsFeatureNamesModForOTU)
e0b5980139d9 maaslin
george-weingart
parents:
diff changeset
168 csvCG.close()