From 38922fcff5a03abfedffda3e06a45fad2270a044 Mon Sep 17 00:00:00 2001 From: westcott Date: Mon, 2 Mar 2009 16:01:29 +0000 Subject: [PATCH] worked on unifrac.unweighted() --- tree.cpp | 53 ++++-------------------------------- tree.h | 5 ++-- unifracunweightedcommand.cpp | 31 ++++++++++----------- unweighted.cpp | 13 +++++---- 4 files changed, 28 insertions(+), 74 deletions(-) diff --git a/tree.cpp b/tree.cpp index 5029438..642c665 100644 --- a/tree.cpp +++ b/tree.cpp @@ -303,12 +303,9 @@ map Tree::mergeGcounts(int position) { } /**************************************************************************************************/ -void Tree::randomLabels() { +void Tree::randomLabels(vector g) { try { - //set up the groups the user wants to include - setGroups(); - for(int i = 0; i < numLeaves; i++){ int z; //get random index to switch with @@ -318,8 +315,8 @@ void Tree::randomLabels() { //if either of the leaf nodes you are about to switch are not in the users groups then you don't want to switch them. bool treez, treei; - treez = inUsersGroups(tree[z].getGroup(), globaldata->Groups); - treei = inUsersGroups(tree[i].getGroup(), globaldata->Groups); + treez = inUsersGroups(tree[z].getGroup(), g); + treei = inUsersGroups(tree[i].getGroup(), g); if ((treez == true) && (treei == true)) { //switches node i and node z's info. @@ -405,8 +402,8 @@ void Tree::randomBlengths() { } } /*************************************************************************************************/ -void Tree::assembleRandomUnifracTree() { - randomLabels(); +void Tree::assembleRandomUnifracTree(vector g) { + randomLabels(g); assembleTree(); } /*************************************************************************************************/ @@ -538,46 +535,6 @@ void Tree::printBranch(int node) { /*****************************************************************/ -void Tree::setGroups() { - try { - //if the user has not entered specific groups to analyze then do them all - if (globaldata->Groups.size() != 0) { - //check that groups are valid - for (int i = 0; i < globaldata->Groups.size(); i++) { - if (globaldata->gTreemap->isValidGroup(globaldata->Groups[i]) != true) { - cout << globaldata->Groups[i] << " is not a valid group, and will be disregarded." << endl; - // erase the invalid group from globaldata->Groups - globaldata->Groups.erase (globaldata->Groups.begin()+i); - } - } - - //if the user only entered invalid groups - if (globaldata->Groups.size() == 0) { - cout << "When using the groups parameter you must have at least 1 valid group. I will run the command using all the groups in your groupfile." << endl; - for (int i = 0; i < globaldata->gTreemap->namesOfGroups.size(); i++) { - globaldata->Groups.push_back(globaldata->gTreemap->namesOfGroups[i]); - } - } - - }else { - for (int i = 0; i < globaldata->gTreemap->namesOfGroups.size(); i++) { - globaldata->Groups.push_back(globaldata->gTreemap->namesOfGroups[i]); - } - } - } - catch(exception& e) { - cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function setGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n"; - exit(1); - } - catch(...) { - cout << "An unknown error has occurred in the Tree class function setGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n"; - exit(1); - } - -} - -/*****************************************************************/ - void Tree::printTree() { for(int i=0;i); void assembleRandomUnifracTree(string, string); void createNewickFile(string); int getIndex(string); @@ -51,11 +51,10 @@ private: map Tree::mergeGcounts(int); void randomTopology(); void randomBlengths(); - void randomLabels(); + void randomLabels(vector); void randomLabels(string, string); int findRoot(); //return index of root node void printBranch(int); //recursively print out tree - void setGroups(); }; #endif diff --git a/unifracunweightedcommand.cpp b/unifracunweightedcommand.cpp index c345ca1..f1326a8 100644 --- a/unifracunweightedcommand.cpp +++ b/unifracunweightedcommand.cpp @@ -85,24 +85,21 @@ int UnifracUnweightedCommand::execute() { //get unweighted scores for random trees for (int j = 0; j < iters; j++) { - int count = 0; - for (int r=0; rgetValues(randT, "", ""); - - //add trees unweighted score to map of scores - it2 = rscoreFreq[count].find(randomData[count]); - if (it2 != rscoreFreq[count].end()) {//already have that score - rscoreFreq[count][randomData[count]]++; - }else{//first time we have seen this score - rscoreFreq[count][randomData[count]] = 1; - } + //we need a different getValues because when we swap the labels we only want to swap those in each parwise comparison + randomData = unweighted->getValues(randT, "", ""); - //add randoms score to validscores - validScores[count][randomData[count]] = randomData[count]; - count++; + for(int k = 0; k < numComp; k++) { +cout << "iter " << j << " comp " << k << " = " << randomData[k] << endl; + //add trees unweighted score to map of scores + it2 = rscoreFreq[k].find(randomData[k]); + if (it2 != rscoreFreq[k].end()) {//already have that score + rscoreFreq[k][randomData[k]]++; + }else{//first time we have seen this score + rscoreFreq[k][randomData[k]] = 1; } + + //add randoms score to validscores + validScores[k][randomData[k]] = randomData[k]; } } @@ -235,7 +232,7 @@ void UnifracUnweightedCommand::setGroups() { } }else { for (int i = 0; i < globaldata->Groups.size(); i++) { - allGroups += tmap->namesOfGroups[i]; + allGroups += globaldata->Groups[i]; numGroups++; } } diff --git a/unweighted.cpp b/unweighted.cpp index b038f38..16a2203 100644 --- a/unweighted.cpp +++ b/unweighted.cpp @@ -205,7 +205,7 @@ EstOutput Unweighted::getValues(Tree* t, string groupA, string groupB) { //if the users enters no groups then give them the score of all groups int numGroups = globaldata->Groups.size(); - +cout << "NumGroups = " << numGroups << endl; //calculate number of comparsions int numComp = 0; for (int r=0; rgetCopy(t); - - //swap labels in the groups you want to compare - copyTree->assembleRandomUnifracTree(globaldata->Groups[a], globaldata->Groups[l]); - + //groups in this combo groups.push_back(globaldata->Groups[a]); groups.push_back(globaldata->Groups[l]); + + //swap labels in the groups you want to compare + copyTree->assembleRandomUnifracTree(groups[0], groups[1]); + for(int i=t->getNumLeaves();igetNumNodes();i++){ @@ -312,7 +313,7 @@ EstOutput Unweighted::getValues(Tree* t, string groupA, string groupB) { copyTree->getCopy(t); //swap labels in all the groups you want to compare - copyTree->assembleRandomUnifracTree(); + copyTree->assembleRandomUnifracTree(groups); for(int i=t->getNumLeaves();igetNumNodes();i++){ -- 2.39.2