]> git.donarmstrong.com Git - mothur.git/commitdiff
fixed bug in bootstrap command that was caused by globaldata's breakup, and made...
authorwestcott <westcott>
Mon, 29 Jun 2009 14:43:21 +0000 (14:43 +0000)
committerwestcott <westcott>
Mon, 29 Jun 2009 14:43:21 +0000 (14:43 +0000)
bootstrapsharedcommand.cpp
bootstrapsharedcommand.h
commandfactory.cpp
concensuscommand.cpp
concensuscommand.h
engine.cpp
mothur.cpp
mothur.h
treegroupscommand.cpp

index ca21a57c4c34e5c09d8de39e8525e972a16d32f7..1fa968f5930243f0efc4d6b0f0764a6892fd1203 100644 (file)
@@ -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<SharedRAbundVector*> 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");
index 4d2ebe39fc3a3283ca2aa6ac844d6e2f5865f9e2..ee393f9eb3fecd8cac0ddfbea3c5e69863b8a621 100644 (file)
@@ -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<Tree*> > trees;  //a vector of trees for each calculator chosen
        vector<Calculator*> treeCalculators;
        vector<ofstream*> out;
        vector< vector<float> > simMatrix;
index db11f8c0dc890dff846787dee0fc48ea8930384a..3ad7425071a938025921cc3f143ccb69781a24ce 100644 (file)
@@ -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);                  }
index 30da32cfe05db5de942cebd213cf9cc9631cb717..37609c3fb4858ca32c02be71539da5547a4c909c 100644 (file)
 
 //**********************************************************************************************************************
 
-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);
index 970102ac27f7a15e6dbfe7b0db2370129861aff3..7b4aa27350857d524181dc554ce1d22f96217b9b 100644 (file)
@@ -38,7 +38,7 @@ private:
        map< vector<string>, int > nodePairsInTree;
        map<string, int>::iterator it;
        map< vector<string>, 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
                                                                                
index 55617bc3f80a97e1832ff61c3e09678de8eca1ca..42938d8a0bf157e0d31a837134d477bc65dca9c3 100644 (file)
@@ -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) {
index 88caa25fd2a16375e6ad515e989cdf19ddcfc916..0cd4faf1dcc3eec95d9cc05b39e32980fe88d5c4 100644 (file)
@@ -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());
index 20d3b4b64725d0428248b47b22a3146a0a3a52c9..4d22e98e2a2bb3eac4dcbbe98b81694b58a00bcb 100644 (file)
--- 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 {
index cb797d8d229fd0a2dfaff82d431b084bd0a227a5..dae04f556490260903f99b4f349d6faea249e26e 100644 (file)
@@ -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();