annotate cpo_galaxy_tree.py @ 6:cabceaa239e4 draft

planemo upload
author jjjjia
date Thu, 23 Aug 2018 12:21:15 -0400
parents fea89c4d5227
children 4d2777aa99db
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
1 #!/home/jjjjia/.conda/envs/py36/bin/python
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
2
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
3 #$ -S /home/jjjjia/.conda/envs/py36/bin/python
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
4 #$ -V # Pass environment variables to the job
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
5 #$ -N CPO_pipeline # Replace with a more specific job name
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
6 #$ -wd /home/jjjjia/testCases # Use the current working dir
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
7 #$ -pe smp 8 # Parallel Environment (how many cores)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
8 #$ -l h_vmem=11G # Memory (RAM) allocation *per core*
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
9 #$ -e ./logs/$JOB_ID.err
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
10 #$ -o ./logs/$JOB_ID.log
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
11 #$ -m ea
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
12 #$ -M bja20@sfu.ca
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
13
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
14 #~/scripts/pipeline.py -i BC11-Kpn005_S2 -f /data/jjjjia/R1/BC11-Kpn005_S2_L001_R1_001.fastq.gz -r /data/jjjjia/R2/BC11-Kpn005_S2_L001_R2_001.fastq.gz -o pipelineResultsQsub -e "Klebsiella pneumoniae"
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
15
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
16 import subprocess
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
17 import pandas
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
18 import optparse
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
19 import os
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
20 import datetime
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
21 import sys
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
22 import time
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
23 import urllib.request
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
24 import gzip
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
25 import collections
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
26 import json
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
27 import numpy
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
28 import ete3 as e
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
29
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
30
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
31 #parses some parameters
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
32 parser = optparse.OptionParser("Usage: %prog [options] arg1 arg2 ...")
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
33 parser.add_option("-t", "--tree", dest="treePath", type="string", default="./pipelineTest/tree.txt", help="identifier of the isolate")
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
34 parser.add_option("-d", "--distance", dest="distancePath", type="string", default="./pipelineTest/distance.tab", help="absolute file path forward read (R1)")
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
35 parser.add_option("-m", "--metadata", dest="metadataPath", type="string", default="./pipelineTest/metadata.tsv",help="absolute file path to reverse read (R2)")
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
36 (options,args) = parser.parse_args()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
37 treePath = str(options.treePath).lstrip().rstrip()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
38 distancePath = str(options.distancePath).lstrip().rstrip()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
39 metadataPath = str(options.metadataPath).lstrip().rstrip()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
40
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
41
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
42 #region result objects
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
43 #define some objects to store values from results
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
44 #//TODO this is not the proper way of get/set private object variables. every value has manually assigned defaults intead of specified in init(). Also, use property(def getVar, def setVar).
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
45 class workflowResult(object):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
46 def __init__(self):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
47 self.new = False
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
48 self.ID = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
49 self.ExpectedSpecies = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
50 self.MLSTSpecies = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
51 self.SequenceType = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
52 self.MLSTScheme = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
53 self.CarbapenemResistanceGenes =""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
54 self.OtherAMRGenes=""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
55 self.TotalPlasmids = 0
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
56 self.plasmids = []
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
57 self.DefinitelyPlasmidContigs =""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
58 self.LikelyPlasmidContigs=""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
59 self.row = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
60 class plasmidObj(object):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
61 def __init__(self):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
62 self.PlasmidsID = 0
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
63 self.Num_Contigs = 0
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
64 self.PlasmidLength = 0
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
65 self.PlasmidRepType = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
66 self.PlasmidMobility = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
67 self.NearestReference = ""
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
68
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
69 #endregion
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
70
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
71 #region useful functions
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
72 def read(path):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
73 return [line.rstrip('\n') for line in open(path)]
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
74 def execute(command):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
75 process = subprocess.Popen(command, shell=False, cwd=curDir, universal_newlines=True, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
76
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
77 # Poll process for new output until finished
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
78 while True:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
79 nextline = process.stdout.readline()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
80 if nextline == '' and process.poll() is not None:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
81 break
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
82 sys.stdout.write(nextline)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
83 sys.stdout.flush()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
84
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
85 output = process.communicate()[0]
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
86 exitCode = process.returncode
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
87
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
88 if (exitCode == 0):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
89 return output
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
90 else:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
91 raise subprocess.CalledProcessError(exitCode, command)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
92 def httpGetFile(url, filepath=""):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
93 if (filepath == ""):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
94 return urllib.request.urlretrieve(url)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
95 else:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
96 urllib.request.urlretrieve(url, filepath)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
97 return True
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
98 def gunzip(inputpath="", outputpath=""):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
99 if (outputpath == ""):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
100 with gzip.open(inputpath, 'rb') as f:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
101 gzContent = f.read()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
102 return gzContent
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
103 else:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
104 with gzip.open(inputpath, 'rb') as f:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
105 gzContent = f.read()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
106 with open(outputpath, 'wb') as out:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
107 out.write(gzContent)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
108 return True
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
109 def addFace(name):
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
110 #if its the reference branch, populate the faces with column headers
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
111 face = e.faces.TextFace(name,fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
112 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
113 face.margin_right = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
114 face.margin_left = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
115 return face
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
116 #endregion
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
117
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
118 #region functions to parse result files
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
119 def ParseWorkflowResults(pathToResult):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
120 _worflowResult = {}
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
121 r = pandas.read_csv(pathToResult, delimiter='\t', header=0)
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
122 r = r.replace(numpy.nan, '', regex=True)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
123 for i in range(len(r.index)):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
124 _results = workflowResult()
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
125 if(str(r.loc[r.index[i], 'new']).lower() == "new"):
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
126 _results.new = True
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
127 else:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
128 _results.new = False
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
129 _results.ID = str(r.loc[r.index[i], 'ID'])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
130 _results.ExpectedSpecies = str(r.loc[r.index[i], 'Expected Species'])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
131 _results.MLSTSpecies = str(r.loc[r.index[i], 'MLST Species'])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
132 _results.SequenceType = str(r.loc[r.index[i], 'Sequence Type'])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
133 _results.MLSTScheme = (str(r.loc[r.index[i], 'MLST Scheme']))
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
134 _results.CarbapenemResistanceGenes = (str(r.loc[r.index[i], 'Carbapenem Resistance Genes']))
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
135 _results.OtherAMRGenes = (str(r.loc[r.index[i], 'Other AMR Genes']))
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
136 _results.TotalPlasmids = int(r.loc[r.index[i], 'Total Plasmids'])
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
137 for j in range(0,_results.TotalPlasmids):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
138 _plasmid = plasmidObj()
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
139 _plasmid.PlasmidsID =(((str(r.loc[r.index[i], 'Plasmids ID'])).split(";"))[j])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
140 _plasmid.Num_Contigs = (((str(r.loc[r.index[i], 'Num_Contigs'])).split(";"))[j])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
141 _plasmid.PlasmidLength = (((str(r.loc[r.index[i], 'Plasmid Length'])).split(";"))[j])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
142 _plasmid.PlasmidRepType = (((str(r.loc[r.index[i], 'Plasmid RepType'])).split(";"))[j])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
143 _plasmid.PlasmidMobility = ((str(r.loc[r.index[i], 'Plasmid Mobility'])).split(";"))[j]
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
144 _plasmid.NearestReference = ((str(r.loc[r.index[i], 'Nearest Reference'])).split(";"))[j]
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
145 _results.plasmids.append(_plasmid)
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
146 _results.DefinitelyPlasmidContigs = (str(r.loc[r.index[i], 'Definitely Plasmid Contigs']))
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
147 _results.LikelyPlasmidContigs = (str(r.loc[r.index[i], 'Likely Plasmid Contigs']))
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
148 _results.row = "\t".join(str(x) for x in r.ix[i].tolist())
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
149 _worflowResult[_results.ID] = _results
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
150 return _worflowResult
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
151
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
152 #endregion
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
153
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
154 def Main():
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
155 metadata = ParseWorkflowResults(metadataPath)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
156 distance = read(distancePath)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
157 treeFile = "".join(read(treePath))
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
158
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
159 distanceDict = {} #store the distance matrix as rowname:list<string>
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
160 for i in range(len(distance)):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
161 temp = distance[i].split("\t")
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
162 distanceDict[temp[0]] = temp[1:]
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
163 #region step5: tree construction
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
164
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
165 '''
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
166 #region create detailed tree
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
167
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
168 plasmidCount = 0
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
169 for n in t.traverse():
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
170 if (n.is_leaf() and not n.name == "Reference"):
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
171 mData = metadata[n.name.replace(".fa","")]
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
172 face = faces.TextFace(mData.MLSTSpecies,fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
173 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
174 face.margin_left = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
175 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
176 n.add_face(face, 0, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
177 face = faces.TextFace(mData.SequenceType,fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
178 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
179 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
180 n.add_face(face, 1, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
181 face = faces.TextFace(mData.CarbapenemResistanceGenes,fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
182 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
183 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
184 n.add_face(face, 2, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
185 index = 3
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
186 if (mData.TotalPlasmids > plasmidCount):
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
187 plasmidCount = mData.TotalPlasmids
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
188 for i in range(0, mData.TotalPlasmids):
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
189 face = faces.TextFace(mData.plasmids[i].PlasmidRepType,fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
190 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
191 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
192 n.add_face(face, index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
193 index+=1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
194 face = faces.TextFace(mData.plasmids[i].PlasmidMobility,fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
195 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
196 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
197 n.add_face(face, index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
198 index+=1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
199
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
200 face = faces.TextFace("Species",fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
201 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
202 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
203 face.margin_left = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
204 (t&"Reference").add_face(face, 0, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
205 face = faces.TextFace("Sequence Type",fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
206 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
207 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
208 (t&"Reference").add_face(face, 1, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
209 face = faces.TextFace("Carbapenamases",fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
210 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
211 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
212 (t&"Reference").add_face(face, 2, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
213 index = 3
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
214 for i in range(0, plasmidCount):
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
215 face = faces.TextFace("plasmid " + str(i) + " replicons",fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
216 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
217 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
218 (t&"Reference").add_face(face, index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
219 index+=1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
220 face = faces.TextFace("plasmid " + str(i) + " mobility",fsize=10,tight_text=True)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
221 face.border.margin = 5
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
222 face.margin_right = 10
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
223 (t&"Reference").add_face(face, index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
224 index+=1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
225
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
226 t.render("./pipelineTest/tree.png", w=5000,units="mm", tree_style=ts)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
227
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
228 #endregion
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
229 '''
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
230 #region create box tree
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
231 #region step5: tree construction
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
232 treeFile = "".join(read(treePath))
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
233 t = e.Tree(treeFile)
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
234 t.set_outgroup(t&"Reference")
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
235
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
236 #set the tree style
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
237 ts = e.TreeStyle()
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
238 ts.show_leaf_name = False
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
239 ts.show_branch_length = True
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
240 ts.scale = 2000 #pixel per branch length unit
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
241 ts.branch_vertical_margin = 15 #pixel between branches
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
242 style2 = e.NodeStyle()
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
243 style2["fgcolor"] = "#000000"
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
244 style2["shape"] = "circle"
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
245 style2["vt_line_color"] = "#0000aa"
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
246 style2["hz_line_color"] = "#0000aa"
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
247 style2["vt_line_width"] = 2
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
248 style2["hz_line_width"] = 2
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
249 style2["vt_line_type"] = 0 # 0 solid, 1 dashed, 2 dotted
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
250 style2["hz_line_type"] = 0
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
251 for n in t.traverse():
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
252 n.set_style(style2)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
253
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
254 #find the plasmid origins
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
255 plasmidIncs = {}
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
256 for key in metadata:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
257 for plasmid in metadata[key].plasmids:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
258 for inc in plasmid.PlasmidRepType.split(","):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
259 if (inc.lower().find("inc") > -1):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
260 if not (inc in plasmidIncs):
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
261 plasmidIncs[inc] = [metadata[key].ID]
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
262 else:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
263 if metadata[key].ID not in plasmidIncs[inc]:
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
264 plasmidIncs[inc].append(metadata[key].ID)
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
265 #plasmidIncs = sorted(plasmidIncs)
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
266 for n in t.traverse(): #loop through the nodes of a tree
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
267 if (n.is_leaf() and n.name == "Reference"):
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
268 #if its the reference branch, populate the faces with column headers
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
269 index = 0
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
270 (t&"Reference").add_face(addFace("SampleID"), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
271 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
272 (t&"Reference").add_face(addFace("New?"), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
273 index = index + 1
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
274 for i in range(len(plasmidIncs)): #this loop adds the columns (aka the incs) to the reference node
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
275 (t&"Reference").add_face(addFace(list(plasmidIncs.keys())[i]), i + index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
276 index = index + len(plasmidIncs)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
277 (t&"Reference").add_face(addFace("MLSTScheme"), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
278 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
279 (t&"Reference").add_face(addFace("Sequence Type"), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
280 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
281 (t&"Reference").add_face(addFace("Carbapenamases"), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
282 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
283 for i in range(len(distanceDict[list(distanceDict.keys())[0]])): #this loop adds the distance matrix
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
284 (t&"Reference").add_face(addFace(distanceDict[list(distanceDict.keys())[0]][i]), index + i, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
285 index = index + len(distanceDict[list(distanceDict.keys())[0]])
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
286 elif (n.is_leaf() and not n.name == "Reference"):
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
287 #not reference branches, populate with metadata
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
288 index = 0
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
289 mData = metadata[n.name.replace(".fa","")]
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
290 n.add_face(addFace(mData.ID), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
291 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
292 if (metadata[n.name.replace(".fa","")].new == True): #new column
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
293 face = e.RectFace(30,30,"green","green") # TextFace("Y",fsize=10,tight_text=True)
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
294 face.border.margin = 5
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
295 face.margin_right = 5
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
296 face.margin_left = 5
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
297 face.vt_align = 1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
298 face.ht_align = 1
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
299 n.add_face(face, index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
300 index = index + 1
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
301 for incs in plasmidIncs: #this loop adds presence/absence to the sample nodes
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
302 if (n.name.replace(".fa","") in plasmidIncs[incs]):
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
303 face = e.RectFace(30,30,"black","black") # TextFace("Y",fsize=10,tight_text=True)
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
304 face.border.margin = 5
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
305 face.margin_right = 5
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
306 face.margin_left = 5
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
307 face.vt_align = 1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
308 face.ht_align = 1
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
309 n.add_face(face, list(plasmidIncs.keys()).index(incs) + index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
310 index = index + len(plasmidIncs)
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
311 n.add_face(addFace(mData.MLSTSpecies), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
312 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
313 n.add_face(addFace(mData.SequenceType), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
314 index = index + 1
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
315 n.add_face(addFace(mData.CarbapenemResistanceGenes), index, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
316 index = index + 1
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
317 for i in range(len(distanceDict[list(distanceDict.keys())[0]])): #this loop adds distance matrix
6
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
318 n.add_face(addFace(list(distanceDict[n.name])[i]), index + i, "aligned")
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
319
cabceaa239e4 planemo upload
jjjjia
parents: 1
diff changeset
320 t.render("./tree.png", w=5000,units="mm", tree_style=ts) #save it as a png. or an phyloxml
1
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
321
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
322 #endregion
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
323 #endregion
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
324
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
325
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
326 start = time.time()#time the analysis
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
327
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
328 #analysis time
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
329 Main()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
330
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
331 end = time.time()
fea89c4d5227 Uploaded
jjjjia
parents:
diff changeset
332 print("Finished!\nThe analysis used: " + str(end-start) + " seconds")