From: westcott Date: Mon, 29 Jun 2009 14:43:21 +0000 (+0000) Subject: fixed bug in bootstrap command that was caused by globaldata's breakup, and made... X-Git-Url: https://git.donarmstrong.com/?p=mothur.git;a=commitdiff_plain;h=0bb069a445785a88a95e3f1427afba7df056a0b3 fixed bug in bootstrap command that was caused by globaldata's breakup, and made concensus command part of bootstrap.shared --- diff --git a/bootstrapsharedcommand.cpp b/bootstrapsharedcommand.cpp index ca21a57..1fa968f 100644 --- a/bootstrapsharedcommand.cpp +++ b/bootstrapsharedcommand.cpp @@ -101,6 +101,9 @@ BootSharedCommand::BootSharedCommand(string option){ if (abort == false) { + //used in tree constructor + globaldata->runParse = false; + validCalculator = new ValidCalculators(); int i; @@ -136,7 +139,10 @@ BootSharedCommand::BootSharedCommand(string option){ for (int i=0; i < treeCalculators.size(); i++) { tempo = new ofstream; out.push_back(tempo); - } + } + + //make a vector of tree* for each calculator + trees.resize(treeCalculators.size()); } } @@ -273,7 +279,9 @@ int BootSharedCommand::execute(){ delete order; } - + + + //reset groups parameter globaldata->Groups.clear(); @@ -286,16 +294,14 @@ int BootSharedCommand::execute(){ } //********************************************************************************************************************** -void BootSharedCommand::createTree(ostream* out){ +void BootSharedCommand::createTree(ostream* out, Tree* t){ try { - //create tree - t = new Tree(); //do merges and create tree structure by setting parents and children //there are numGroups - 1 merges to do for (int i = 0; i < (numGroups - 1); i++) { - float largest = -1.0; + float largest = -1000.0; int row, column; //find largest value in sims matrix by searching lower triangle for (int j = 1; j < simMatrix.size(); j++) { @@ -328,8 +334,8 @@ void BootSharedCommand::createTree(ostream* out){ index[column] = numGroups+i; //zero out highest value that caused the merge. - simMatrix[row][column] = -1.0; - simMatrix[column][row] = -1.0; + simMatrix[row][column] = -1000.0; + simMatrix[column][row] = -1000.0; //merge values in simsMatrix for (int n = 0; n < simMatrix.size(); n++) { @@ -337,10 +343,14 @@ void BootSharedCommand::createTree(ostream* out){ simMatrix[row][n] = (simMatrix[row][n] + simMatrix[column][n]) / 2; simMatrix[n][row] = simMatrix[row][n]; //delete column - simMatrix[column][n] = -1.0; - simMatrix[n][column] = -1.0; + simMatrix[column][n] = -1000.0; + simMatrix[n][column] = -1000.0; } } + + //adjust tree to make sure root to tip length is .5 + int root = t->findRoot(); + t->tree[root].setBranchLength((0.5 - t->tree[root].getLengthToLeaves())); //assemble tree t->assembleTree(); @@ -348,9 +358,6 @@ void BootSharedCommand::createTree(ostream* out){ //print newick file t->print(*out); - //delete tree - delete t; - } catch(exception& e) { errorOut(e, "BootSharedCommand", "createTree"); @@ -380,6 +387,7 @@ void BootSharedCommand::process(SharedOrderVector* order) { EstOutput data; vector subset; + //open an ostream for each calc to print to for (int z = 0; z < treeCalculators.size(); z++) { //create a new filename @@ -387,10 +395,13 @@ void BootSharedCommand::process(SharedOrderVector* order) { openOutputFile(outputFile, *(out[z])); } + mothurOut("Generating bootstrap trees..."); cout.flush(); + //create a file for each calculator with the 1000 trees in it. for (int p = 0; p < iters; p++) { - util->getSharedVectorswithReplacement(Groups, lookup, order); //fills group vectors from order vector. + util->getSharedVectorswithReplacement(globaldata->Groups, lookup, order); //fills group vectors from order vector. + //for each calculator for(int i = 0 ; i < treeCalculators.size(); i++) { @@ -423,14 +434,47 @@ void BootSharedCommand::process(SharedOrderVector* order) { } } } - + + tempTree = new Tree(); + //creates tree from similarity matrix and write out file - createTree(out[i]); + createTree(out[i], tempTree); + + //save trees for concensus command. + trees[i].push_back(tempTree); } } + + mothurOut("\tDone."); mothurOutEndLine(); + //delete globaldata's tree + for (int m = 0; m < globaldata->gTree.size(); m++) { delete globaldata->gTree[m]; } + globaldata->clear(); + + + //create concensus trees for each bootstrapped tree set + for (int k = 0; k < trees.size(); k++) { + + mothurOut("Generating concensus tree for " + treeCalculators[k]->getName()); mothurOutEndLine(); + + //set global data to calc trees + globaldata->gTree = trees[k]; + + string filename = getRootName(globaldata->inputFileName) + treeCalculators[k]->getName() + ".boot" + order->getLabel(); + concensus = new ConcensusCommand(filename); + concensus->execute(); + delete concensus; + + //delete globaldata's tree + for (int m = 0; m < globaldata->gTree.size(); m++) { delete globaldata->gTree[m]; } + globaldata->clear(); + + } + + + //close ostream for each calc for (int z = 0; z < treeCalculators.size(); z++) { out[z]->close(); } - + } catch(exception& e) { errorOut(e, "BootSharedCommand", "process"); diff --git a/bootstrapsharedcommand.h b/bootstrapsharedcommand.h index 4d2ebe3..ee393f9 100644 --- a/bootstrapsharedcommand.h +++ b/bootstrapsharedcommand.h @@ -19,6 +19,7 @@ #include "tree.h" #include "treemap.h" #include "sharedutilities.h" +#include "concensuscommand.h" class GlobalData; @@ -31,7 +32,7 @@ public: void help(); private: - void createTree(ostream*); + void createTree(ostream*, Tree*); void printSims(); void process(SharedOrderVector*); @@ -41,6 +42,9 @@ private: ReadOTUFile* read; TreeMap* tmap; Tree* t; + Tree* tempTree; + ConcensusCommand* concensus; + vector< vector > trees; //a vector of trees for each calculator chosen vector treeCalculators; vector out; vector< vector > simMatrix; diff --git a/commandfactory.cpp b/commandfactory.cpp index db11f8c..3ad7425 100644 --- a/commandfactory.cpp +++ b/commandfactory.cpp @@ -12,7 +12,6 @@ #include "readtreecommand.h" #include "readotucommand.h" #include "clustercommand.h" -#include "parselistcommand.h" #include "collectcommand.h" #include "collectsharedcommand.h" #include "getgroupcommand.h" @@ -39,7 +38,7 @@ #include "getoturepcommand.h" #include "treegroupscommand.h" #include "bootstrapsharedcommand.h" -#include "concensuscommand.h" +//#include "concensuscommand.h" #include "distancecommand.h" #include "aligncommand.h" #include "matrixoutputcommand.h" @@ -89,7 +88,7 @@ CommandFactory::CommandFactory(){ commands["get.sabund"] = "get.sabund"; commands["get.rabund"] = "get.rabund"; commands["bootstrap.shared"] = "bootstrap.shared"; - commands["concensus"] = "concensus"; + //commands["concensus"] = "concensus"; commands["help"] = "help"; commands["filter.seqs"] = "filter.seqs"; commands["align.seqs"] = "align.seqs"; @@ -146,7 +145,7 @@ Command* CommandFactory::getCommand(string commandName, string optionString){ else if(commandName == "tree.shared") { command = new TreeGroupCommand(optionString); } else if(commandName == "dist.shared") { command = new MatrixOutputCommand(optionString); } else if(commandName == "bootstrap.shared") { command = new BootSharedCommand(optionString); } - else if(commandName == "concensus") { command = new ConcensusCommand(optionString); } + //else if(commandName == "concensus") { command = new ConcensusCommand(optionString); } else if(commandName == "dist.seqs") { command = new DistanceCommand(optionString); } else if(commandName == "align.seqs") { command = new AlignCommand(optionString); } else if(commandName == "summary.seqs") { command = new SeqSummaryCommand(optionString); } diff --git a/concensuscommand.cpp b/concensuscommand.cpp index 30da32c..37609c3 100644 --- a/concensuscommand.cpp +++ b/concensuscommand.cpp @@ -11,21 +11,24 @@ //********************************************************************************************************************** -ConcensusCommand::ConcensusCommand(string option){ +ConcensusCommand::ConcensusCommand(string fileroot){ try { globaldata = GlobalData::getInstance(); abort = false; + filename = fileroot; //allow user to run help - if(option == "help") { help(); abort = true; } + //if(option == "help") { help(); abort = true; } - else { - if (option != "") { mothurOut("There are no valid parameters for the concensus command."); mothurOutEndLine(); abort = true; } + //else { + //if (option != "") { mothurOut("There are no valid parameters for the concensus command."); mothurOutEndLine(); abort = true; } - //no trees were read - if (globaldata->gTree.size() == 0) { mothurOut("You must execute the read.tree command, before you may use the concensus command."); mothurOutEndLine(); abort = true; } - else { t = globaldata->gTree; } - } + // //no trees were read + // if (globaldata->gTree.size() == 0) { mothurOut("You must execute the read.tree command, before you may use the concensus command."); mothurOutEndLine(); abort = true; } + //else { + t = globaldata->gTree; + // } + //} } catch(exception& e) { errorOut(e, "ConcensusCommand", "ConcensusCommand"); @@ -71,7 +74,7 @@ int ConcensusCommand::execute(){ getSets(); //open file for pairing not included in the tree - notIncluded = getRootName(globaldata->inputFileName) + "concensuspairs"; + notIncluded = filename + ".concensuspairs"; openOutputFile(notIncluded, out2); concensusTree = new Tree(); @@ -140,7 +143,7 @@ int ConcensusCommand::execute(){ out2 << '\t' << it2->second << endl; } - outputFile = getRootName(globaldata->inputFileName) + "concensus.tre"; + outputFile = filename + ".cons.tre"; openOutputFile(outputFile, out); concensusTree->printForBoot(out); diff --git a/concensuscommand.h b/concensuscommand.h index 970102a..7b4aa27 100644 --- a/concensuscommand.h +++ b/concensuscommand.h @@ -38,7 +38,7 @@ private: map< vector, int > nodePairsInTree; map::iterator it; map< vector, int>::iterator it2; - string outputFile, notIncluded; + string outputFile, notIncluded, filename; ofstream out, out2; int numNodes, numLeaves, count; //count is the next available spot in the tree vector diff --git a/engine.cpp b/engine.cpp index 55617bc..42938d8 100644 --- a/engine.cpp +++ b/engine.cpp @@ -21,7 +21,7 @@ InteractEngine::InteractEngine(string path){ globaldata = GlobalData::getInstance(); globaldata->argv = path; - system("clear"); + } /***********************************************************************/ @@ -86,7 +86,7 @@ BatchEngine::BatchEngine(string path, string batchFileName){ openedBatch = openInputFile(batchFileName, inputBatchFile); globaldata->argv = path; - system("clear"); + // char buffer = ' '; // ifstream header("introtext.txt"); @@ -182,7 +182,7 @@ ScriptEngine::ScriptEngine(string path, string commandString){ globaldata->argv = path; - system("clear"); + } catch(exception& e) { diff --git a/mothur.cpp b/mothur.cpp index 88caa25..0cd4faf 100644 --- a/mothur.cpp +++ b/mothur.cpp @@ -18,6 +18,8 @@ GlobalData* GlobalData::_uniqueInstance = 0; int main(int argc, char *argv[]){ try { + system("clear"); + //remove old logfile string logFileName = "mothur.logFile"; remove(logFileName.c_str()); diff --git a/mothur.h b/mothur.h index 20d3b4b..4d22e98 100644 --- a/mothur.h +++ b/mothur.h @@ -153,14 +153,13 @@ string toString(const T&x, int i){ return output.str(); } - /***********************************************************************/ inline int openOutputFileAppend(string fileName, ofstream& fileHandle){ fileHandle.open(fileName.c_str(), ios::app); if(!fileHandle) { - cerr << "Error: Could not open " << fileName << endl; + cout << "Error: Could not open " << fileName << endl; return 1; } else { @@ -169,7 +168,6 @@ inline int openOutputFileAppend(string fileName, ofstream& fileHandle){ } - /**************************************************************************************************/ inline void mothurOut(string message) { @@ -239,6 +237,8 @@ inline void errorOut(exception& e, string object, string function) { } + + /***********************************************************************/ inline void gobble(istream& f){ @@ -394,7 +394,7 @@ inline int openInputFile(string fileName, ifstream& fileHandle){ fileHandle.open(fileName.c_str()); if(!fileHandle) { - cerr << "Error: Could not open " << fileName << endl; + mothurOut("Error: Could not open " + fileName); mothurOutEndLine(); return 1; } else { @@ -409,7 +409,7 @@ inline int openOutputFile(string fileName, ofstream& fileHandle){ fileHandle.open(fileName.c_str(), ios::trunc); if(!fileHandle) { - cerr << "Error: Could not open " << fileName << endl; + mothurOut("Error: Could not open " + fileName); mothurOutEndLine(); return 1; } else { diff --git a/treegroupscommand.cpp b/treegroupscommand.cpp index cb797d8..dae04f5 100644 --- a/treegroupscommand.cpp +++ b/treegroupscommand.cpp @@ -227,7 +227,8 @@ int TreeGroupCommand::execute(){ lastLabel = lookup[0]->getLabel(); if (lookup.size() < 2) { mothurOut("You have not provided enough valid groups. I cannot run the command."); mothurOutEndLine(); return 0; } - + + //used in tree constructor globaldata->runParse = false; //create tree file @@ -266,6 +267,7 @@ int TreeGroupCommand::execute(){ //fills globaldatas tree names globaldata->Treenames = globaldata->Groups; + //used in tree constructor globaldata->runParse = false; makeSimsDist();