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 Weighted : public TreeCalculator {
22 Weighted( bool r) : includeRoot(r) {};
25 EstOutput getValues(Tree*, string, string);
26 EstOutput getValues(Tree*, int, string);
32 linePair(int i, int j) : start(i), num(j) {}
34 vector<linePair> lines;
37 map<string, int>::iterator it;
38 map<string, double> WScore; //a score for each group combination i.e. AB, AC, BC.
41 map< vector<string>, set<int> > rootForGrouping; //maps a grouping combo to the root for that combo
44 EstOutput driver(Tree*, vector< vector<string> >, int, int, CountTable*);
45 EstOutput createProcesses(Tree*, vector< vector<string> >, CountTable*);
46 double getLengthToRoot(Tree*, int, string, string);
49 /***********************************************************************/
55 vector< vector<string> > namesOfGroupCombos;
62 weightedData(MothurOut* mout, int st, int en, vector< vector<string> > ngc, Tree* tree, CountTable* count, bool ir) {
66 namesOfGroupCombos = ngc;
73 /**************************************************************************************************/
74 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux) || (__linux__) || (__unix__) || (__unix)
76 static DWORD WINAPI MyWeightedThreadFunction(LPVOID lpParam){
77 weightedData* pDataArray;
78 pDataArray = (weightedData*)lpParam;
80 map<string, int>::iterator it;
83 map< vector<string>, set<int> > rootForGrouping;
84 map<string, double> WScore;
86 for (int h = pDataArray->start; h < (pDataArray->start+pDataArray->num); h++) {
88 if (pDataArray->m->control_pressed) { return 0; }
90 //initialize weighted score
91 string groupA = pDataArray->namesOfGroupCombos[h][0];
92 string groupB = pDataArray->namesOfGroupCombos[h][1];
94 set<int> validBranches;
95 WScore[groupA+groupB] = 0.0;
96 D.push_back(0.0000); //initialize a spot in D for each combination
98 //adding the wieghted sums from group i
99 for (int j = 0; j < pDataArray->t->groupNodeInfo[groupA].size(); j++) { //the leaf nodes that have seqs from group i
100 map<string, int>::iterator it = pDataArray->t->tree[pDataArray->t->groupNodeInfo[groupA][j]].pcount.find(groupA);
101 int numSeqsInGroupI = it->second;
103 //double sum = getLengthToRoot(pDataArray->t, pDataArray->t->groupNodeInfo[groupA][j], groupA, groupB);
104 /*************************************************************************************/
106 int index = pDataArray->t->groupNodeInfo[groupA][j];
109 if(pDataArray->t->tree[index].getBranchLength() != -1){ sum += abs(pDataArray->t->tree[index].getBranchLength()); }
110 double tempTotal = 0.0;
111 index = pDataArray->t->tree[index].getParent();
113 vector<string> grouping; grouping.push_back(groupA); grouping.push_back(groupB);
115 rootForGrouping[grouping].insert(index);
117 //while you aren't at root
118 while(pDataArray->t->tree[index].getParent() != -1){
120 if (pDataArray->m->control_pressed) { return 0; }
122 int parent = pDataArray->t->tree[index].getParent();
124 if (pDataArray->includeRoot) { //add everyone
125 if(pDataArray->t->tree[index].getBranchLength() != -1){ sum += abs(pDataArray->t->tree[index].getBranchLength()); }
128 //am I the root for this grouping? if so I want to stop "early"
129 //does my sibling have descendants from the users groups?
130 int lc = pDataArray->t->tree[parent].getLChild();
131 int rc = pDataArray->t->tree[parent].getRChild();
134 if (lc == index) { sib = rc; }
136 map<string, int>::iterator itGroup;
138 itGroup = pDataArray->t->tree[sib].pcount.find(groupA);
139 if (itGroup != pDataArray->t->tree[sib].pcount.end()) { pcountSize++; }
140 itGroup = pDataArray->t->tree[sib].pcount.find(groupB);
141 if (itGroup != pDataArray->t->tree[sib].pcount.end()) { pcountSize++; }
143 //if yes, I am not the root so add me
144 if (pcountSize != 0) {
145 if (pDataArray->t->tree[index].getBranchLength() != -1) {
146 sum += abs(pDataArray->t->tree[index].getBranchLength()) + tempTotal;
152 rootForGrouping[grouping].clear();
153 rootForGrouping[grouping].insert(parent);
154 }else { //if no, I may be the root so add my br to tempTotal until I am proven innocent
155 if (pDataArray->t->tree[index].getBranchLength() != -1) {
156 tempTotal += abs(pDataArray->t->tree[index].getBranchLength());
164 //get all nodes above the root to add so we don't add their u values above
165 index = *(rootForGrouping[grouping].begin());
167 while(pDataArray->t->tree[index].getParent() != -1){
168 int parent = pDataArray->t->tree[index].getParent();
169 rootForGrouping[grouping].insert(parent);
173 /*************************************************************************************/
174 double weightedSum = ((numSeqsInGroupI * sum) / (double)pDataArray->ct->getGroupCount(groupA));
176 D[count] += weightedSum;
179 //adding the wieghted sums from group l
180 for (int j = 0; j < pDataArray->t->groupNodeInfo[groupB].size(); j++) { //the leaf nodes that have seqs from group l
181 map<string, int>::iterator it = pDataArray->t->tree[pDataArray->t->groupNodeInfo[groupB][j]].pcount.find(groupB);
182 int numSeqsInGroupL = it->second;
184 //double sum = getLengthToRoot(pDataArray->t, pDataArray->t->groupNodeInfo[groupB][j], groupA, groupB);
185 /*************************************************************************************/
187 int index = pDataArray->t->groupNodeInfo[groupB][j];
190 if(pDataArray->t->tree[index].getBranchLength() != -1){ sum += abs(pDataArray->t->tree[index].getBranchLength()); }
191 double tempTotal = 0.0;
192 index = pDataArray->t->tree[index].getParent();
194 vector<string> grouping; grouping.push_back(groupA); grouping.push_back(groupB);
196 rootForGrouping[grouping].insert(index);
198 //while you aren't at root
199 while(pDataArray->t->tree[index].getParent() != -1){
201 if (pDataArray->m->control_pressed) { return 0; }
203 int parent = pDataArray->t->tree[index].getParent();
205 if (pDataArray->includeRoot) { //add everyone
206 if(pDataArray->t->tree[index].getBranchLength() != -1){ sum += abs(pDataArray->t->tree[index].getBranchLength()); }
209 //am I the root for this grouping? if so I want to stop "early"
210 //does my sibling have descendants from the users groups?
211 int lc = pDataArray->t->tree[parent].getLChild();
212 int rc = pDataArray->t->tree[parent].getRChild();
215 if (lc == index) { sib = rc; }
217 map<string, int>::iterator itGroup;
219 itGroup = pDataArray->t->tree[sib].pcount.find(groupA);
220 if (itGroup != pDataArray->t->tree[sib].pcount.end()) { pcountSize++; }
221 itGroup = pDataArray->t->tree[sib].pcount.find(groupB);
222 if (itGroup != pDataArray->t->tree[sib].pcount.end()) { pcountSize++; }
224 //if yes, I am not the root so add me
225 if (pcountSize != 0) {
226 if (pDataArray->t->tree[index].getBranchLength() != -1) {
227 sum += abs(pDataArray->t->tree[index].getBranchLength()) + tempTotal;
233 rootForGrouping[grouping].clear();
234 rootForGrouping[grouping].insert(parent);
235 }else { //if no, I may be the root so add my br to tempTotal until I am proven innocent
236 if (pDataArray->t->tree[index].getBranchLength() != -1) {
237 tempTotal += abs(pDataArray->t->tree[index].getBranchLength());
245 //get all nodes above the root to add so we don't add their u values above
246 index = *(rootForGrouping[grouping].begin());
248 while(pDataArray->t->tree[index].getParent() != -1){
249 int parent = pDataArray->t->tree[index].getParent();
250 rootForGrouping[grouping].insert(parent);
254 /*************************************************************************************/
256 double weightedSum = ((numSeqsInGroupL * sum) / (double)pDataArray->ct->getGroupCount(groupB));
258 D[count] += weightedSum;
263 //calculate u for the group comb
264 for (int h = pDataArray->start; h < (pDataArray->start+pDataArray->num); h++) {
266 //pDataArray->m->mothurOut("Processing combo: " + toString(h)); pDataArray->m->mothurOutEndLine();
268 string groupA = pDataArray->namesOfGroupCombos[h][0];
269 string groupB = pDataArray->namesOfGroupCombos[h][1];
271 //calculate u for the group comb
272 for(int i=0;i<pDataArray->t->getNumNodes();i++){
274 if (pDataArray->m->control_pressed) { return 0; }
277 //int pcountSize = 0;
278 //does this node have descendants from groupA
279 it = pDataArray->t->tree[i].pcount.find(groupA);
280 //if it does u = # of its descendants with a certain group / total number in tree with a certain group
281 if (it != pDataArray->t->tree[i].pcount.end()) {
282 u = (double) pDataArray->t->tree[i].pcount[groupA] / (double) pDataArray->ct->getGroupCount(groupA);
286 //does this node have descendants from group l
287 it = pDataArray->t->tree[i].pcount.find(groupB);
289 //if it does subtract their percentage from u
290 if (it != pDataArray->t->tree[i].pcount.end()) {
291 u -= (double) pDataArray->t->tree[i].pcount[groupB] / (double) pDataArray->ct->getGroupCount(groupB);
294 if (pDataArray->includeRoot) {
295 if (pDataArray->t->tree[i].getBranchLength() != -1) {
296 u = abs(u * pDataArray->t->tree[i].getBranchLength());
297 WScore[(groupA+groupB)] += u;
300 //if this is not the root then add it
301 if (rootForGrouping[pDataArray->namesOfGroupCombos[h]].count(i) == 0) {
302 if (pDataArray->t->tree[i].getBranchLength() != -1) {
303 u = abs(u * pDataArray->t->tree[i].getBranchLength());
304 WScore[(groupA+groupB)] += u;
312 /********************************************************/
313 //calculate weighted score for the group combination
316 for (int h = pDataArray->start; h < (pDataArray->start+pDataArray->num); h++) {
317 UN = (WScore[pDataArray->namesOfGroupCombos[h][0]+pDataArray->namesOfGroupCombos[h][1]] / D[count]);
318 if (isnan(UN) || isinf(UN)) { UN = 0; }
319 pDataArray->results.push_back(UN);
326 catch(exception& e) {
327 pDataArray->m->errorOut(e, "Weighted", "MyWeightedThreadFunction");