weighted.cpp

   1 /*
   2  *  weighted.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 "weighted.h"
  11
  12 /**************************************************************************************************/
  13
  14 EstOutput Weighted::getValues(Tree* t) {
  15     try {
  16
  17                 int numGroups = tmap->getNumGroups();
  18
  19                 //calculate number of comparisons i.e. with groups A,B,C = AB, AC, BC = 3;
  20                 int n = 1;
  21                 for (int i=1; i<numGroups; i++) {
  22                         for (int l = n; l < numGroups; l++) {
  23                                 //initialize weighted scores
  24                                 WScore[tmap->namesOfGroups[i-1]+tmap->namesOfGroups[l]] = 0.0;
  25                         }
  26                 }
  27
  28                 data.clear(); //clear out old values
  29
  30                 double D = 0.0000;
  31
  32                 for(int i=0;i<t->getNumLeaves();i++){
  33                         int index = i;
  34                         double sum = 0.0000;
  35
  36                         //while you aren't at root
  37                         while(t->tree[index].getParent() != -1){
  38
  39                                 if(t->tree[index].pGroups.size() != 0){
  40                                         sum += t->tree[index].getBranchLength();
  41                                 }
  42
  43                                 //old_index = you
  44                                 int old_index = index;
  45                                 //index = your parent
  46                                 index = t->tree[index].getParent();
  47
  48                                 //if you grandparent is the root
  49                                 if(t->tree[index].getParent() == -1 && t->tree[old_index].pGroups.size() != 0){
  50                                         int lc = t->tree[t->tree[index].getLChild()].pGroups.size();
  51                                         int rc = t->tree[t->tree[index].getRChild()].pGroups.size();
  52
  53
  54                                         if(lc == 0 || rc == 0){
  55                                                 sum -= t->tree[old_index].getBranchLength();
  56                                         }
  57                                 }
  58                         }
  59
  60                         if(t->tree[i].getGroup() != ""){
  61                                 sum /= (double)tmap->seqsPerGroup[t->tree[i].getGroup()];
  62                                 D += sum;
  63                         }
  64                 }
  65
  66
  67                 for(int i=0;i<t->getNumNodes();i++){
  68                         //calculate weighted score for each of the group comb i.e. with groups A,B,C = AB, AC, BC.
  69                         n = 1;
  70                         for (int b=1; b<numGroups; b++) {
  71                                 for (int l = n; l < numGroups; l++) {
  72                                         double u;
  73                                         //does this node have descendants from group b-1
  74                                         it = t->tree[i].pcount.find(tmap->namesOfGroups[b-1]);
  75                                         //if it does u = # of its descendants with a certain group / total number in tree with a certain group
  76                                         if (it != t->tree[i].pcount.end()) {
  77                                                 u = (double) t->tree[i].pcount[tmap->namesOfGroups[b-1]] / (double) tmap->seqsPerGroup[tmap->namesOfGroups[b-1]];
  78                                         }else { u = 0.00; }
  79
  80                                         //does this node have descendants from group l
  81                                         it = t->tree[i].pcount.find(tmap->namesOfGroups[l]);
  82                                         //if it does subtract their percentage from u
  83                                         if (it != t->tree[i].pcount.end()) {
  84                                                 u -= (double) t->tree[i].pcount[tmap->namesOfGroups[l]] / (double) tmap->seqsPerGroup[tmap->namesOfGroups[l]];
  85                                         }
  86
  87                                         u = abs(u) * t->tree[i].getBranchLength();
  88
  89                                         //save groupcombs u value
  90                                         WScore[tmap->namesOfGroups[b-1]+tmap->namesOfGroups[l]] += u;
  91
  92                                 }
  93                                 n++;
  94                         }
  95                 }
  96
  97                 //calculate weighted score for each group combination
  98                 double UN;
  99                 n = 1;
 100                 for (int i=1; i<numGroups; i++) {
 101                         for (int l = n; l < numGroups; l++) {
 102                                 UN = (WScore[tmap->namesOfGroups[i-1]+tmap->namesOfGroups[l]] / D);
 103                                 if (isnan(UN) || isinf(UN)) { UN = 0; }
 104                                 data.push_back(UN);
 105                         }
 106                 }
 107
 108                 return data;
 109         }
 110         catch(exception& e) {
 111                 cout << "Standard Error: " << e.what() << " has occurred in the Weighted class Function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 112                 exit(1);
 113         }
 114         catch(...) {
 115                 cout << "An unknown error has occurred in the Weighted class function getValues. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 116                 exit(1);
 117         }
 118
 119 }
 120
 121 /**************************************************************************************************/