9 * Created by Sarah Westcott on 2/9/09.
10 * Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
14 #include "treecalculator.h"
15 #include "counttable.h"
17 /***********************************************************************/
19 class Unweighted : public TreeCalculator {
22 Unweighted(bool r) : includeRoot(r) {};
24 EstOutput getValues(Tree*, int, string);
25 EstOutput getValues(Tree*, string, string, int, string);
31 linePair(int i, int j) : start(i), num(j) {}
33 vector<linePair> lines;
38 map< vector<string>, set<int> > rootForGrouping; //maps a grouping combo to the roots for that combo
41 EstOutput driver(Tree*, vector< vector<string> >, int, int, CountTable*);
42 EstOutput createProcesses(Tree*, vector< vector<string> >, CountTable*);
43 EstOutput driver(Tree*, vector< vector<string> >, int, int, bool, CountTable*);
44 EstOutput createProcesses(Tree*, vector< vector<string> >, bool, CountTable*);
45 int getRoot(Tree*, int, vector<string>);
48 /***********************************************************************/
49 struct unweightedData {
54 vector< vector<string> > namesOfGroupCombos;
60 unweightedData(MothurOut* mout, int st, int en, vector< vector<string> > ngc, Tree* tree, CountTable* count, bool ir) {
64 namesOfGroupCombos = ngc;
71 /**************************************************************************************************/
72 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
74 static DWORD WINAPI MyUnWeightedThreadFunction(LPVOID lpParam){
75 unweightedData* pDataArray;
76 pDataArray = (unweightedData*)lpParam;
78 pDataArray->results.resize(pDataArray->num);
79 map< vector<string>, set<int> > rootForGrouping;
83 for (int h = pDataArray->start; h < (pDataArray->start+pDataArray->num); h++) {
84 if (pDataArray->m->control_pressed) { return 0; }
86 double UniqueBL=0.0000; //a branch length is unique if it's chidren are from the same group
87 double totalBL = 0.00; //all branch lengths
88 double UW = 0.00; //Unweighted Value = UniqueBL / totalBL;
90 //find a node that belongs to one of the groups in this combo
91 int nodeBelonging = -1;
92 for (int g = 0; g < pDataArray->namesOfGroupCombos[h].size(); g++) {
93 if (pDataArray->t->groupNodeInfo[pDataArray->namesOfGroupCombos[h][g]].size() != 0) { nodeBelonging = pDataArray->t->groupNodeInfo[pDataArray->namesOfGroupCombos[h][g]][0]; break; }
97 if (nodeBelonging == -1) {
98 pDataArray->m->mothurOut("[WARNING]: cannot find a nodes in the tree from grouping ");
99 for (int g = 0; g < pDataArray->namesOfGroupCombos[h].size()-1; g++) { pDataArray->m->mothurOut(pDataArray->namesOfGroupCombos[h][g] + "-"); }
100 pDataArray->m->mothurOut(pDataArray->namesOfGroupCombos[h][pDataArray->namesOfGroupCombos[h].size()-1]);
101 pDataArray->m->mothurOut(", skipping."); pDataArray->m->mothurOutEndLine(); pDataArray->results[count] = UW;
104 //if including the root this clears rootForGrouping[namesOfGroupCombos[h]]
105 //getRoot(t, nodeBelonging, namesOfGroupCombos[h]);
106 /////////////////////////////////////////////////////////////////////////////
107 //you are a leaf so get your parent
108 vector<string> grouping = pDataArray->namesOfGroupCombos[h];
109 int index = pDataArray->t->tree[nodeBelonging].getParent();
111 if (pDataArray->includeRoot) {
112 rootForGrouping[grouping].clear();
115 //my parent is a potential root
116 rootForGrouping[grouping].insert(index);
118 //while you aren't at root
119 while(pDataArray->t->tree[index].getParent() != -1){
120 //cout << index << endl;
121 if (pDataArray->m->control_pressed) { return 0; }
123 //am I the root for this grouping? if so I want to stop "early"
124 //does my sibling have descendants from the users groups?
125 //if so I am not the root
126 int parent = pDataArray->t->tree[index].getParent();
127 int lc = pDataArray->t->tree[parent].getLChild();
128 int rc = pDataArray->t->tree[parent].getRChild();
131 if (lc == index) { sib = rc; }
133 map<string, int>::iterator itGroup;
135 for (int j = 0; j < grouping.size(); j++) {
136 map<string, int>::iterator itGroup = pDataArray->t->tree[sib].pcount.find(grouping[j]);
137 if (itGroup != pDataArray->t->tree[sib].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
140 //if yes, I am not the root
141 if (pcountSize != 0) {
142 rootForGrouping[grouping].clear();
143 rootForGrouping[grouping].insert(parent);
149 //get all nodes above the root to add so we don't add their u values above
150 index = *(rootForGrouping[grouping].begin());
151 while(pDataArray->t->tree[index].getParent() != -1){
152 int parent = pDataArray->t->tree[index].getParent();
153 rootForGrouping[grouping].insert(parent);
154 //cout << parent << " in root" << endl;
158 /////////////////////////////////////////////////////////////////////////////
160 for(int i=0;i<pDataArray->t->getNumNodes();i++){
162 if (pDataArray->m->control_pressed) { return 0; }
164 //pcountSize = 0, they are from a branch that is entirely from a group the user doesn't want
165 //pcountSize = 2, not unique to one group
166 //pcountSize = 1, unique to one group
169 for (int j = 0; j < pDataArray->namesOfGroupCombos[h].size(); j++) {
170 map<string, int>::iterator itGroup = pDataArray->t->tree[i].pcount.find(pDataArray->namesOfGroupCombos[h][j]);
171 if (itGroup != pDataArray->t->tree[i].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
176 if (pcountSize == 0) { }
177 else if ((pDataArray->t->tree[i].getBranchLength() != -1) && (pcountSize == 1) && (rootForGrouping[pDataArray->namesOfGroupCombos[h]].count(i) == 0)) { //you have a unique branch length and you are not the root
178 UniqueBL += abs(pDataArray->t->tree[i].getBranchLength());
182 if (pcountSize == 0) { }
183 else if ((pDataArray->t->tree[i].getBranchLength() != -1) && (pcountSize != 0) && (rootForGrouping[pDataArray->namesOfGroupCombos[h]].count(i) == 0)) { //you have a branch length and you are not the root
184 totalBL += abs(pDataArray->t->tree[i].getBranchLength());
187 //cout << UniqueBL << '\t' << totalBL << endl;
188 UW = (UniqueBL / totalBL);
190 if (isnan(UW) || isinf(UW)) { UW = 0; }
192 pDataArray->results[count] = UW;
200 catch(exception& e) {
201 pDataArray->m->errorOut(e, "UnWeighted", "MyUnWeightedThreadFunction");
205 /**************************************************************************************************/
207 static DWORD WINAPI MyUnWeightedRandomThreadFunction(LPVOID lpParam){
208 unweightedData* pDataArray;
209 pDataArray = (unweightedData*)lpParam;
211 pDataArray->results.resize(pDataArray->num);
215 Tree* copyTree = new Tree(pDataArray->ct);
217 for (int h = pDataArray->start; h < (pDataArray->start+pDataArray->num); h++) {
219 if (pDataArray->m->control_pressed) { return 0; }
221 map< vector<string>, set<int> > rootForGrouping;
223 //copy random tree passed in
224 copyTree->getCopy(pDataArray->t);
226 //swap labels in the groups you want to compare
227 copyTree->assembleRandomUnifracTree(pDataArray->namesOfGroupCombos[h]);
229 double UniqueBL=0.0000; //a branch length is unique if it's chidren are from the same group
230 double totalBL = 0.00; //all branch lengths
231 double UW = 0.00; //Unweighted Value = UniqueBL / totalBL;
232 //find a node that belongs to one of the groups in this combo
233 int nodeBelonging = -1;
234 for (int g = 0; g < pDataArray->namesOfGroupCombos[h].size(); g++) {
235 if (copyTree->groupNodeInfo[pDataArray->namesOfGroupCombos[h][g]].size() != 0) { nodeBelonging = copyTree->groupNodeInfo[pDataArray->namesOfGroupCombos[h][g]][0]; break; }
239 if (nodeBelonging == -1) {
240 pDataArray->m->mothurOut("[WARNING]: cannot find a nodes in the tree from grouping ");
241 for (int g = 0; g < pDataArray->namesOfGroupCombos[h].size()-1; g++) { pDataArray->m->mothurOut(pDataArray->namesOfGroupCombos[h][g] + "-"); }
242 pDataArray->m->mothurOut(pDataArray->namesOfGroupCombos[h][pDataArray->namesOfGroupCombos[h].size()-1]);
243 pDataArray->m->mothurOut(", skipping."); pDataArray->m->mothurOutEndLine(); pDataArray->results[count] = UW;
246 //if including the root this clears rootForGrouping[namesOfGroupCombos[h]]
247 //getRoot(copyTree, nodeBelonging, namesOfGroupCombos[h]);
248 /////////////////////////////////////////////////////////////////////////////
249 //you are a leaf so get your parent
250 vector<string> grouping = pDataArray->namesOfGroupCombos[h];
251 int index = copyTree->tree[nodeBelonging].getParent();
253 if (pDataArray->includeRoot) {
254 rootForGrouping[grouping].clear();
257 //my parent is a potential root
258 rootForGrouping[grouping].insert(index);
260 //while you aren't at root
261 while(copyTree->tree[index].getParent() != -1){
262 //cout << index << endl;
263 if (pDataArray->m->control_pressed) { return 0; }
265 //am I the root for this grouping? if so I want to stop "early"
266 //does my sibling have descendants from the users groups?
267 //if so I am not the root
268 int parent = copyTree->tree[index].getParent();
269 int lc = copyTree->tree[parent].getLChild();
270 int rc = copyTree->tree[parent].getRChild();
273 if (lc == index) { sib = rc; }
275 map<string, int>::iterator itGroup;
277 for (int j = 0; j < grouping.size(); j++) {
278 map<string, int>::iterator itGroup = copyTree->tree[sib].pcount.find(grouping[j]);
279 if (itGroup != copyTree->tree[sib].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
282 //if yes, I am not the root
283 if (pcountSize != 0) {
284 rootForGrouping[grouping].clear();
285 rootForGrouping[grouping].insert(parent);
291 //get all nodes above the root to add so we don't add their u values above
292 index = *(rootForGrouping[grouping].begin());
293 while(copyTree->tree[index].getParent() != -1){
294 int parent = copyTree->tree[index].getParent();
295 rootForGrouping[grouping].insert(parent);
296 //cout << parent << " in root" << endl;
300 /////////////////////////////////////////////////////////////////////////////
301 for(int i=0;i<copyTree->getNumNodes();i++){
303 if (pDataArray->m->control_pressed) { return 0; }
305 //pcountSize = 0, they are from a branch that is entirely from a group the user doesn't want
306 //pcountSize = 2, not unique to one group
307 //pcountSize = 1, unique to one group
309 for (int j = 0; j < pDataArray->namesOfGroupCombos[h].size(); j++) {
310 map<string, int>::iterator itGroup = copyTree->tree[i].pcount.find(pDataArray->namesOfGroupCombos[h][j]);
311 if (itGroup != copyTree->tree[i].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
315 if (pcountSize == 0) { }
316 else if ((copyTree->tree[i].getBranchLength() != -1) && (pcountSize == 1) && (rootForGrouping[pDataArray->namesOfGroupCombos[h]].count(i) == 0)) { //you have a unique branch length and you are not the root
317 UniqueBL += abs(copyTree->tree[i].getBranchLength());
321 if (pcountSize == 0) { }
322 else if ((copyTree->tree[i].getBranchLength() != -1) && (pcountSize != 0) && (rootForGrouping[pDataArray->namesOfGroupCombos[h]].count(i) == 0)) { //you have a branch length and you are not the root
323 totalBL += abs(copyTree->tree[i].getBranchLength());
327 cout << h << '\t' << UniqueBL << '\t' << totalBL << endl;
328 UW = (UniqueBL / totalBL);
330 if (isnan(UW) || isinf(UW)) { UW = 0; }
332 pDataArray->results[count] = UW;
333 cout << h << '\t' << UW << endl;
343 catch(exception& e) {
344 pDataArray->m->errorOut(e, "UnWeighted", "MyUnWeightedRandomThreadFunction");