]> git.donarmstrong.com Git - mothur.git/blob - unweighted.cpp
worked on unifrac.unweighted()
[mothur.git] / unweighted.cpp
1 /*
2  *  unweighted.cpp
3  *  Mothur
4  *
5  *  Created by Sarah Westcott on 2/9/09.
6  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
7  *
8  */
9
10 #include "unweighted.h"
11
12 /**************************************************************************************************/
13
14 EstOutput Unweighted::getValues(Tree* t) {
15         try {
16                 globaldata = GlobalData::getInstance();
17                 
18                 vector<string> groups;
19                 double UniqueBL;  //a branch length is unique if it's chidren are from the same group
20                 double totalBL; //all branch lengths
21                 double UW;              //Unweighted Value = UniqueBL / totalBL;
22                 map<string, int>::iterator it;  //iterator to traverse pgroups
23                 map<string, int> copyIpcount;
24
25         
26                 //if the users enters no groups then give them the score of all groups
27                 int numGroups = globaldata->Groups.size();
28                 
29                 //calculate number of comparsions
30                 int numComp = 0;
31                 for (int r=0; r<numGroups; r++) { 
32                         for (int l = r+1; l < numGroups; l++) {
33                                 numComp++;
34                         }
35                 }
36
37                 //numComp+1 for AB, AC, BC, ABC
38                 data.resize(numComp+1,0);
39                 
40                 int count = 0;
41                 for (int a=0; a<numGroups; a++) { 
42                         for (int l = a+1; l < numGroups; l++) {
43                                 UniqueBL=0.0000;  //a branch length is unique if it's chidren are from the same group
44                                 totalBL = 0.00; //all branch lengths
45                                 UW = 0.00;              //Unweighted Value = UniqueBL / totalBL;
46                                 copyIpcount.clear();
47                                 
48                                 //groups in this combo
49                                 groups.push_back(globaldata->Groups[a]); groups.push_back(globaldata->Groups[l]);
50                 
51                                 for(int i=t->getNumLeaves();i<t->getNumNodes();i++){
52                 
53                                         int lc = t->tree[i].getLChild();  //lc = vector index of left child
54                                         int rc = t->tree[i].getRChild();  //rc = vector index of right child
55                         
56                                         /**********************************************************************/
57                                         //This section adds in all lengths that are non leaf
58                         
59                                         copyIpcount = t->tree[i].pcount;
60                                         for (it = copyIpcount.begin(); it != copyIpcount.end(); it++) {
61                                                 if (inUsersGroups(it->first, groups) != true) { copyIpcount.erase(it->first);   }
62                                         }
63                         
64                                         //if i's children are from the same group then i's pcount size will be 1 
65                                         //if copyIpcount.size() = 0 they are from a branch that is entirely from a group the user doesn't want
66                                         if (copyIpcount.size() == 0) { }
67                                         else if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() == 1)) {  UniqueBL += t->tree[i].getBranchLength();        }
68                         
69                                         //add i's BL to total if it is from the groups the user wants
70                                         if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() != 0)) {  
71                                                 totalBL += t->tree[i].getBranchLength(); 
72                                         }
73                         
74                                         /**********************************************************************/
75                                         //This section adds in all lengths that are leaf
76                         
77                                         //if i's chidren are leaves
78                                         if (t->tree[rc].getRChild() == -1) {
79                                                 //if rc is a valid group and rc has a BL
80                                                 if ((inUsersGroups(t->tree[rc].getGroup(), groups) == true) && (t->tree[rc].getBranchLength() != -1)) {
81                                                         UniqueBL += t->tree[rc].getBranchLength();
82                                                         totalBL += t->tree[rc].getBranchLength(); 
83                                                 }
84                                         }
85                         
86                                         if (t->tree[lc].getLChild() == -1) {
87                                                 //if lc is a valid group and lc has a BL
88                                                 if ((inUsersGroups(t->tree[lc].getGroup(), groups) == true) && (t->tree[lc].getBranchLength() != -1)) {
89                                                         UniqueBL += t->tree[lc].getBranchLength();
90                                                         totalBL += t->tree[lc].getBranchLength(); 
91                                                 }
92                                         }
93                         
94                                         /**********************************************************************/
95                                 }
96                 
97                                 UW = (UniqueBL / totalBL);  
98         
99                                 if (isnan(UW) || isinf(UW)) { UW = 0; }
100         
101                                 data[count] = UW;
102                                 count++;
103                                 groups.clear();
104                         }
105                 }
106                 
107                 
108                 if (numComp != 1) {
109                         if (numGroups == 0) {
110                                 //get score for all users groups
111                                 for (int i = 0; i < tmap->namesOfGroups.size(); i++) {
112                                         groups.push_back(tmap->namesOfGroups[i]);
113                                 }
114                         }else {
115                                 for (int i = 0; i < globaldata->Groups.size(); i++) {
116                                         groups.push_back(globaldata->Groups[i]);
117                                 }
118                         }
119                 
120                         UniqueBL=0.0000;  //a branch length is unique if it's chidren are from the same group
121                         totalBL = 0.00; //all branch lengths
122                         UW = 0.00;              //Unweighted Value = UniqueBL / totalBL;
123                         copyIpcount.clear();
124                                 
125                         for(int i=t->getNumLeaves();i<t->getNumNodes();i++){
126                 
127                                 int lc = t->tree[i].getLChild();  //lc = vector index of left child
128                                 int rc = t->tree[i].getRChild();  //rc = vector index of right child
129                         
130                                 /**********************************************************************/
131                                 //This section adds in all lengths that are non leaf
132                         
133                                 copyIpcount = t->tree[i].pcount;
134                                 for (it = copyIpcount.begin(); it != copyIpcount.end(); it++) {
135                                         if (inUsersGroups(it->first, groups) != true) { copyIpcount.erase(it->first);   }
136                                 }
137                         
138                                 //if i's children are from the same group then i's pcount size will be 1 
139                                 //if copyIpcount.size() = 0 they are from a branch that is entirely from a group the user doesn't want
140                                 if (copyIpcount.size() == 0) { }
141                                 else if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() == 1)) {  UniqueBL += t->tree[i].getBranchLength();        }
142                         
143                                 //add i's BL to total if it is from the groups the user wants
144                                 if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() != 0)) {  
145                                         totalBL += t->tree[i].getBranchLength(); 
146                                 }
147                         
148                                 /**********************************************************************/
149                                 //This section adds in all lengths that are leaf
150                         
151                                 //if i's chidren are leaves
152                                 if (t->tree[rc].getRChild() == -1) {
153                                         //if rc is a valid group and rc has a BL
154                                         if ((inUsersGroups(t->tree[rc].getGroup(), groups) == true) && (t->tree[rc].getBranchLength() != -1)) {
155                                                 UniqueBL += t->tree[rc].getBranchLength();
156                                                 totalBL += t->tree[rc].getBranchLength(); 
157                                         }
158                                 }
159                         
160                                 if (t->tree[lc].getLChild() == -1) {
161                                         //if lc is a valid group and lc has a BL
162                                         if ((inUsersGroups(t->tree[lc].getGroup(), groups) == true) && (t->tree[lc].getBranchLength() != -1)) {
163                                                 UniqueBL += t->tree[lc].getBranchLength();
164                                                 totalBL += t->tree[lc].getBranchLength(); 
165                                         }
166                                 }
167                         
168                                 /**********************************************************************/
169                         }
170                 
171                         UW = (UniqueBL / totalBL);  
172         
173                         if (isnan(UW) || isinf(UW)) { UW = 0; }
174         
175                         data[count] = UW;
176                 }
177
178                 return data;
179         
180         }
181         catch(exception& e) {
182                 cout << "Standard Error: " << e.what() << " has occurred in the Unweighted class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
183                 exit(1);
184         }
185         catch(...) {
186                 cout << "An unknown error has occurred in the Unweighted class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
187                 exit(1);
188         }
189
190 }
191
192 /**************************************************************************************************/
193
194 EstOutput Unweighted::getValues(Tree* t, string groupA, string groupB) { 
195  try {
196         globaldata = GlobalData::getInstance();
197                 
198                 vector<string> groups;
199                 double UniqueBL;  //a branch length is unique if it's chidren are from the same group
200                 double totalBL; //all branch lengths
201                 double UW;              //Unweighted Value = UniqueBL / totalBL;
202                 map<string, int>::iterator it;  //iterator to traverse pgroups
203                 map<string, int> copyIpcount;
204                 copyTree = new Tree;
205
206                 //if the users enters no groups then give them the score of all groups
207                 int numGroups = globaldata->Groups.size();
208 cout << "NumGroups = " << numGroups << endl;
209                 //calculate number of comparsions
210                 int numComp = 0;
211                 for (int r=0; r<numGroups; r++) { 
212                         for (int l = r+1; l < numGroups; l++) {
213                                 numComp++;
214                         }
215                 }
216
217                 //numComp+1 for AB, AC, BC, ABC
218                 data.resize(numComp+1,0);
219                 
220                 int count = 0;
221                 for (int a=0; a<numGroups; a++) { 
222                         for (int l = a+1; l < numGroups; l++) {
223                                 UniqueBL=0.0000;  //a branch length is unique if it's chidren are from the same group
224                                 totalBL = 0.00; //all branch lengths
225                                 UW = 0.00;              //Unweighted Value = UniqueBL / totalBL;
226                                 copyIpcount.clear();
227                                 
228                                 //copy random tree passed in
229                                 copyTree->getCopy(t);
230                                                                 
231                                 //groups in this combo
232                                 groups.push_back(globaldata->Groups[a]); groups.push_back(globaldata->Groups[l]);
233                                 
234                                 //swap labels in the groups you want to compare
235                                 copyTree->assembleRandomUnifracTree(groups[0], groups[1]);
236
237                 
238                                 for(int i=t->getNumLeaves();i<t->getNumNodes();i++){
239                 
240                                         int lc = t->tree[i].getLChild();  //lc = vector index of left child
241                                         int rc = t->tree[i].getRChild();  //rc = vector index of right child
242                         
243                                         /**********************************************************************/
244                                         //This section adds in all lengths that are non leaf
245                         
246                                         copyIpcount = t->tree[i].pcount;
247                                         for (it = copyIpcount.begin(); it != copyIpcount.end(); it++) {
248                                                 if (inUsersGroups(it->first, groups) != true) { copyIpcount.erase(it->first);   }
249                                         }
250                         
251                                         //if i's children are from the same group then i's pcount size will be 1 
252                                         //if copyIpcount.size() = 0 they are from a branch that is entirely from a group the user doesn't want
253                                         if (copyIpcount.size() == 0) { }
254                                         else if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() == 1)) {  UniqueBL += t->tree[i].getBranchLength();        }
255                         
256                                         //add i's BL to total if it is from the groups the user wants
257                                         if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() != 0)) {  
258                                                 totalBL += t->tree[i].getBranchLength(); 
259                                         }
260                         
261                                         /**********************************************************************/
262                                         //This section adds in all lengths that are leaf
263                         
264                                         //if i's chidren are leaves
265                                         if (t->tree[rc].getRChild() == -1) {
266                                                 //if rc is a valid group and rc has a BL
267                                                 if ((inUsersGroups(t->tree[rc].getGroup(), groups) == true) && (t->tree[rc].getBranchLength() != -1)) {
268                                                         UniqueBL += t->tree[rc].getBranchLength();
269                                                         totalBL += t->tree[rc].getBranchLength(); 
270                                                 }
271                                         }
272                         
273                                         if (t->tree[lc].getLChild() == -1) {
274                                                 //if lc is a valid group and lc has a BL
275                                                 if ((inUsersGroups(t->tree[lc].getGroup(), groups) == true) && (t->tree[lc].getBranchLength() != -1)) {
276                                                         UniqueBL += t->tree[lc].getBranchLength();
277                                                         totalBL += t->tree[lc].getBranchLength(); 
278                                                 }
279                                         }
280                         
281                                         /**********************************************************************/
282                                 }
283                 
284                                 UW = (UniqueBL / totalBL);  
285         
286                                 if (isnan(UW) || isinf(UW)) { UW = 0; }
287         
288                                 data[count] = UW;
289                                 count++;
290                                 groups.clear();
291                         }
292                 }
293                 
294                 
295                 if (numComp != 1) {
296                         if (numGroups == 0) {
297                                 //get score for all users groups
298                                 for (int i = 0; i < tmap->namesOfGroups.size(); i++) {
299                                         groups.push_back(tmap->namesOfGroups[i]);
300                                 }
301                         }else {
302                                 for (int i = 0; i < globaldata->Groups.size(); i++) {
303                                         groups.push_back(globaldata->Groups[i]);
304                                 }
305                         }
306                 
307                         UniqueBL=0.0000;  //a branch length is unique if it's chidren are from the same group
308                         totalBL = 0.00; //all branch lengths
309                         UW = 0.00;              //Unweighted Value = UniqueBL / totalBL;
310                         copyIpcount.clear();
311                 
312                         //copy random tree passed in
313                         copyTree->getCopy(t);
314                                 
315                         //swap labels in all the groups you want to compare
316                         copyTree->assembleRandomUnifracTree(groups);
317
318                         for(int i=t->getNumLeaves();i<t->getNumNodes();i++){
319                 
320                                 int lc = t->tree[i].getLChild();  //lc = vector index of left child
321                                 int rc = t->tree[i].getRChild();  //rc = vector index of right child
322                         
323                                 /**********************************************************************/
324                                 //This section adds in all lengths that are non leaf
325                         
326                                 copyIpcount = t->tree[i].pcount;
327                                 for (it = copyIpcount.begin(); it != copyIpcount.end(); it++) {
328                                         if (inUsersGroups(it->first, groups) != true) { copyIpcount.erase(it->first);   }
329                                 }
330                         
331                                 //if i's children are from the same group then i's pcount size will be 1 
332                                 //if copyIpcount.size() = 0 they are from a branch that is entirely from a group the user doesn't want
333                                 if (copyIpcount.size() == 0) { }
334                                 else if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() == 1)) {  UniqueBL += t->tree[i].getBranchLength();        }
335                         
336                                 //add i's BL to total if it is from the groups the user wants
337                                 if ((t->tree[i].getBranchLength() != -1) && (copyIpcount.size() != 0)) {  
338                                         totalBL += t->tree[i].getBranchLength(); 
339                                 }
340                         
341                                 /**********************************************************************/
342                                 //This section adds in all lengths that are leaf
343                         
344                                 //if i's chidren are leaves
345                                 if (t->tree[rc].getRChild() == -1) {
346                                         //if rc is a valid group and rc has a BL
347                                         if ((inUsersGroups(t->tree[rc].getGroup(), groups) == true) && (t->tree[rc].getBranchLength() != -1)) {
348                                                 UniqueBL += t->tree[rc].getBranchLength();
349                                                 totalBL += t->tree[rc].getBranchLength(); 
350                                         }
351                                 }
352                         
353                                 if (t->tree[lc].getLChild() == -1) {
354                                         //if lc is a valid group and lc has a BL
355                                         if ((inUsersGroups(t->tree[lc].getGroup(), groups) == true) && (t->tree[lc].getBranchLength() != -1)) {
356                                                 UniqueBL += t->tree[lc].getBranchLength();
357                                                 totalBL += t->tree[lc].getBranchLength(); 
358                                         }
359                                 }
360                         
361                                 /**********************************************************************/
362                         }
363                 
364                         UW = (UniqueBL / totalBL);  
365         
366                         if (isnan(UW) || isinf(UW)) { UW = 0; }
367         
368                         data[count] = UW;
369                 }
370
371                 return data;
372         
373         }
374         catch(exception& e) {
375                 cout << "Standard Error: " << e.what() << " has occurred in the Unweighted class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
376                 exit(1);
377         }
378         catch(...) {
379                 cout << "An unknown error has occurred in the Unweighted class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
380                 exit(1);
381         }
382 }
383
384
385
386