Entering edit mode
3.3 years ago
1165492577
•
0
I've got error in executing the python script provided by CAFE
[XX@XX DIR]$ cafe cafetutorial_run.sh
# output: resultfile.cafe
(python27) [XX@XX DIR]$ python cafetutorial_report_analysis.py -in resultfile.cafe -o summary_run
=======================================================================
CAFE Report File Analysis
01.22.2021 | 10:14:46
---------
Parsing format information...
---------
Initializing output structures...
---------
Counting changes per branch...
[*-------------------------------------------------] 0.135% complete.Traceback (most recent call last):
File "/new_data/SOFTWARE/phylogenetics/CAFE_python_scripts/cafetutorial_report_analysis.py", line 371, in <module>
results_main, node_fams_main = cra(inlines_main, results_main, node_fams_main, linestart_main, ancfilename, sorted_nodes, 1);
File "/new_data/SOFTWARE/phylogenetics/CAFE_python_scripts/cafetutorial_report_analysis.py", line 192, in cra
curcount = int(tlinfo[tlnode][4]);
ValueError: invalid literal for int() with base 10: '6_0'
Here is part of the python script:
131 #######################
132 def cra(inlines, results, node_fams, linestart, afilename, s_nodes, v):
133 numbars = 0;
134 donepercent = [];
135 i = 0;
136 acount = 0;
137 afile = open(afilename, "a");
138
139 for inline in inlines:
140 # Each line of the report file is read.
141 if v == 1:
142 numbars, donepercent = cafecore.loadingBar(i, len(inlines), donepercent, numbars);
143 i = i + 1;
144
145 if i <= linestart:
146 continue;
147 # If the line is a CAFE info line, skip it.
148
149 inline = inline.strip().split("\t");
150 famid = inline[0];
151 famtree = inline[1];
152 nodeformat = inline[3].replace("),(", ") (");
153 # Parsing the information for the current family.
154
155 outline = famid + "\t";
156 outlist = [0 for n in s_nodes];
157 # Prep for the anc states file.
158
159 nodes = formatLineParse(nodeformat);
160
161 tlinfo, newfamtree = cafecore.treeParseNew(famtree, 1);
162 # Reading the tree and adding my node labels.
163
164 for tlnode in tlinfo:
165 if tlinfo[tlnode][3] == 'root':
166 tlinfo[tlnode].append(famtree[famtree.rfind("_")+1:]);
167 elif tlinfo[tlnode][3] == 'tip':
168 tlinfo[tlnode].append(tlnode[tlnode.index("_")+1:]);
169 else:
170 tlinfo[tlnode][4] = tlinfo[tlnode][4][1:];
171 # Gene counts for each node are read as support values for internal nodes, but must
172 # have the underscore removed. Tip and root node counts are added here as well.
173
174 tlinfo = nodeRelabel(tlinfo);
175 # Removes the gene counts from the tip node labels.
176
177 if i == (linestart + 1):
178 maps = nodeMap(tinfo, tlinfo);
179 # If this is the first family, we need to build our maps from my node ids to CAFE's.
180
181 for tlnode in tlinfo:
182 # For each node in the current gene family tree, we make our counts.
183
184 if tlinfo[tlnode][3] == 'root':
185 continue;
186 # No count is made at the root of the tree.
187
188 curanc = tlinfo[tlnode][1];
189 curmap = maps[tlnode];
190 # Get the ancestor and the map of the current node.
191
192 curcount = int(tlinfo[tlnode][4]);
193 anccount = int(tlinfo[curanc][4]);
194 # Get the gene counts of the current node and the ancestor.
195
196 outlist[s_nodes.index(curmap)] = str(curcount);
197 # Save the count of the current node to be sent to the anc states file.
198
199 diff = curcount - anccount;
200 # Calculate the difference in gene count.
201
202 typeflag = 0;
203 # typeflag tells whether an expansion or contraction has occurred.
204
205 if curcount > anccount:
206 typeflag = 1;
207 results[curmap][0] += 1;
208 results[curmap][1] += diff;
209
210 if r_opt == 0:
211 node_fams[curmap][0].append(famid + "[+" + str(diff) + "]");
212
213 if famid not in node_fams['total']:
214 node_fams['total'].append(famid);
215 # If the difference in gene count between the current node and the ancestor is positive, an
216 # expansion has occurred. This makes the appropriate counts.
217
218 elif curcount < anccount:
219 typeflag = 2
220 results[curmap][3] += 1;
221 results[curmap][4] += abs(diff);
222
223 if curcount == 0 and anccount != 0:
224 results[curmap][5] += 1;
225
226 if r_opt == 0:
227 node_fams[curmap][1].append(famid + "[" + str(diff) + "]");
228
229 if famid not in node_fams['total']:
230 node_fams['total'].append(famid);
231 # If the difference in gene count between the current node and the ancestor is negative, a
232 # contraction has occurred. This makes the appropriate counts. It also checks for family losses
233 # along that branch by seeing if the current node has transitioned to a count of 0.
234
235 elif curcount == anccount:
236 results[curmap][2] += 1;
237 # Otherwise, the counts at the current node and the ancestor are the same and no change has occurred.
