From 38922fcff5a03abfedffda3e06a45fad2270a044 Mon Sep 17 00:00:00 2001
From: westcott <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<string,int> Tree::mergeGcounts(int position) {
}
/**************************************************************************************************/
-void Tree::randomLabels() {
+void Tree::randomLabels(vector<string> 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<string> 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<numNodes;i++){
diff --git a/tree.h b/tree.h
index 733f3a0..d6a3b55 100644
--- a/tree.h
+++ b/tree.h
@@ -25,7 +25,7 @@ public:
void getCopy(Tree*); //makes tree a copy of the one passed in.
void assembleRandomTree();
- void assembleRandomUnifracTree();
+ void assembleRandomUnifracTree(vector<string>);
void assembleRandomUnifracTree(string, string);
void createNewickFile(string);
int getIndex(string);
@@ -51,11 +51,10 @@ private:
map<string,int> Tree::mergeGcounts(int);
void randomTopology();
void randomBlengths();
- void randomLabels();
+ void randomLabels(vector<string>);
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; r<numGroups; r++) {
- for (int l = r+1; l < numGroups; l++) {
- //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 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; r<numGroups; r++) {
@@ -227,12 +227,13 @@ EstOutput Unweighted::getValues(Tree* t, string groupA, string groupB) {
//copy random tree passed in
copyTree->getCopy(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();i<t->getNumNodes();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();i<t->getNumNodes();i++){
--
2.39.5