]> git.donarmstrong.com Git - mothur.git/blobdiff - unifracweightedcommand.cpp
added set.current and get.current commands and modified existing commands to update...
[mothur.git] / unifracweightedcommand.cpp
index df73c9389b99b7ba18725d418d62b36184efa398..d5212009b46003b4046c5f8e2895321edbcb907d 100644 (file)
@@ -12,7 +12,7 @@
 //**********************************************************************************************************************
 vector<string> UnifracWeightedCommand::getValidParameters(){   
        try {
-               string Array[] =  {"groups","iters","distance","random","processors","outputdir","inputdir"};
+               string Array[] =  {"groups","iters","distance","random","processors","root","outputdir","inputdir"};
                vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
                return myArray;
        }
@@ -24,12 +24,12 @@ vector<string> UnifracWeightedCommand::getValidParameters(){
 //**********************************************************************************************************************
 UnifracWeightedCommand::UnifracWeightedCommand(){      
        try {
-               abort = true;
-               //initialize outputTypes
+               abort = true; calledHelp = true; 
                vector<string> tempOutNames;
                outputTypes["weighted"] = tempOutNames;
                outputTypes["wsummary"] = tempOutNames;
                outputTypes["phylip"] = tempOutNames;
+               outputTypes["column"] = tempOutNames;
        }
        catch(exception& e) {
                m->errorOut(e, "UnifracWeightedCommand", "UnifracWeightedCommand");
@@ -64,15 +64,15 @@ vector<string> UnifracWeightedCommand::getRequiredFiles(){
 UnifracWeightedCommand::UnifracWeightedCommand(string option) {
        try {
                globaldata = GlobalData::getInstance();
-               abort = false;
+               abort = false; calledHelp = false;   
                Groups.clear();
                        
                //allow user to run help
-               if(option == "help") { help(); abort = true; }
+               if(option == "help") { help(); abort = true; calledHelp = true; }
                
                else {
                        //valid paramters for this command
-                       string Array[] =  {"groups","iters","distance","random","processors","outputdir","inputdir"};
+                       string Array[] =  {"groups","iters","distance","random","processors","root","outputdir","inputdir"};
                        vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
                        
                        OptionParser parser(option);
@@ -90,6 +90,7 @@ UnifracWeightedCommand::UnifracWeightedCommand(string option) {
                        outputTypes["weighted"] = tempOutNames;
                        outputTypes["wsummary"] = tempOutNames;
                        outputTypes["phylip"] = tempOutNames;
+                       outputTypes["column"] = tempOutNames;
                        
                        if (globaldata->gTree.size() == 0) {//no trees were read
                                m->mothurOut("You must execute the read.tree command, before you may execute the unifrac.weighted command."); m->mothurOutEndLine(); abort = true;  }
@@ -112,12 +113,19 @@ UnifracWeightedCommand::UnifracWeightedCommand(string option) {
                        itersString = validParameter.validFile(parameters, "iters", false);                     if (itersString == "not found") { itersString = "1000"; }
                        convert(itersString, iters); 
                        
-                       string temp = validParameter.validFile(parameters, "distance", false);                  if (temp == "not found") { temp = "false"; }
-                       phylip = m->isTrue(temp);
-               
-                       temp = validParameter.validFile(parameters, "random", false);                                   if (temp == "not found") { temp = "F"; }
+                       string temp = validParameter.validFile(parameters, "distance", false);                  
+                       if (temp == "not found") { phylip = false; outputForm = ""; }
+                       else{
+                               if ((temp == "lt") || (temp == "column") || (temp == "square")) {  phylip = true;  outputForm = temp; }
+                               else { m->mothurOut("Options for distance are: lt, square, or column. Using lt."); m->mothurOutEndLine(); phylip = true; outputForm = "lt"; }
+                       }
+                       
+                       temp = validParameter.validFile(parameters, "random", false);                           if (temp == "not found") { temp = "F"; }
                        random = m->isTrue(temp);
                        
+                       temp = validParameter.validFile(parameters, "root", false);                                     if (temp == "not found") { temp = "F"; }
+                       includeRoot = m->isTrue(temp);
+                       
                        temp = validParameter.validFile(parameters, "processors", false);       if (temp == "not found"){       temp = "1";                             }
                        convert(temp, processors); 
                        
@@ -136,7 +144,7 @@ UnifracWeightedCommand::UnifracWeightedCommand(string option) {
                                util->setGroups(globaldata->Groups, tmap->namesOfGroups, s, numGroups, "weighted");     //sets the groups the user wants to analyze
                                util->getCombos(groupComb, globaldata->Groups, numComp);
                                
-                               weighted = new Weighted(tmap);
+                               weighted = new Weighted(tmap, includeRoot);
                                
                        }
                }
@@ -153,11 +161,12 @@ UnifracWeightedCommand::UnifracWeightedCommand(string option) {
 void UnifracWeightedCommand::help(){
        try {
                m->mothurOut("The unifrac.weighted command can only be executed after a successful read.tree command.\n");
-               m->mothurOut("The unifrac.weighted command parameters are groups, iters, distance, processors and random.  No parameters are required.\n");
+               m->mothurOut("The unifrac.weighted command parameters are groups, iters, distance, processors, root and random.  No parameters are required.\n");
                m->mothurOut("The groups parameter allows you to specify which of the groups in your groupfile you would like analyzed.  You must enter at least 2 valid groups.\n");
                m->mothurOut("The group names are separated by dashes.  The iters parameter allows you to specify how many random trees you would like compared to your tree.\n");
                m->mothurOut("The distance parameter allows you to create a distance file from the results. The default is false.\n");
                m->mothurOut("The random parameter allows you to shut off the comparison to random trees. The default is false, meaning don't compare your trees with randomly generated trees.\n");
+               m->mothurOut("The root parameter allows you to include the entire root in your calculations. The default is false, meaning stop at the root for this comparision instead of the root of the entire tree.\n");
                m->mothurOut("The processors parameter allows you to specify the number of processors to use. The default is 1.\n");
                m->mothurOut("The unifrac.weighted command should be in the following format: unifrac.weighted(groups=yourGroups, iters=yourIters).\n");
                m->mothurOut("Example unifrac.weighted(groups=A-B-C, iters=500).\n");
@@ -175,13 +184,15 @@ void UnifracWeightedCommand::help(){
 int UnifracWeightedCommand::execute() {
        try {
        
-               if (abort == true) { return 0; }
+               if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
                
                int start = time(NULL);
                
                //get weighted for users tree
                userData.resize(numComp,0);  //data[0] = weightedscore AB, data[1] = weightedscore AC...
                randomData.resize(numComp,0); //data[0] = weightedscore AB, data[1] = weightedscore AC...
+               
+               if (numComp < processors) { processors = numComp; }
                                
                //get weighted scores for users trees
                for (int i = 0; i < T.size(); i++) {
@@ -256,7 +267,7 @@ int UnifracWeightedCommand::execute() {
                                        if (m->control_pressed) { delete output; outSum.close(); for (int i = 0; i < outputNames.size(); i++) { remove(outputNames[i].c_str());  } return 0; }
                                        
                                        //report progress
-                                       m->mothurOut("Iter: " + toString(j+1)); m->mothurOutEndLine();          
+//                                     m->mothurOut("Iter: " + toString(j+1)); m->mothurOutEndLine();          
                                }
                                lines.clear();
                        
@@ -307,6 +318,19 @@ int UnifracWeightedCommand::execute() {
                
                m->mothurOut("It took " + toString(time(NULL) - start) + " secs to run unifrac.weighted."); m->mothurOutEndLine();
                
+               //set phylip file as new current phylipfile
+               string current = "";
+               itTypes = outputTypes.find("phylip");
+               if (itTypes != outputTypes.end()) {
+                       if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setPhylipFile(current); }
+               }
+               
+               //set column file as new current columnfile
+               itTypes = outputTypes.find("column");
+               if (itTypes != outputTypes.end()) {
+                       if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setColumnFile(current); }
+               }
+               
                m->mothurOutEndLine();
                m->mothurOut("Output File Names: "); m->mothurOutEndLine();
                for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
@@ -326,7 +350,6 @@ int UnifracWeightedCommand::createProcesses(Tree* t, vector< vector<string> > na
        try {
 #if defined (__APPLE__) || (__MACH__) || (linux) || (__linux)
                int process = 1;
-               int num = 0;
                vector<int> processIDS;
                
                EstOutput results;
@@ -497,15 +520,23 @@ void UnifracWeightedCommand::createPhylipFile() {
                //for each tree
                for (int i = 0; i < T.size(); i++) { 
                
-                       string phylipFileName = outputDir + m->getSimpleName(globaldata->getTreeFile())  + toString(i+1) + ".weighted.dist";
-                       outputNames.push_back(phylipFileName);
-                       outputTypes["phylip"].push_back(phylipFileName);
+                       string phylipFileName;
+                       if ((outputForm == "lt") || (outputForm == "square")) {
+                               phylipFileName = outputDir + m->getSimpleName(globaldata->getTreeFile())  + toString(i+1) + ".weighted.phylip.dist";
+                               outputNames.push_back(phylipFileName); outputTypes["phylip"].push_back(phylipFileName); 
+                       }else { //column
+                               phylipFileName = outputDir + m->getSimpleName(globaldata->getTreeFile())  + toString(i+1) + ".weighted.column.dist";
+                               outputNames.push_back(phylipFileName); outputTypes["column"].push_back(phylipFileName); 
+                       }
+                       
                        ofstream out;
                        m->openOutputFile(phylipFileName, out);
                        
-                       //output numSeqs
-                       out << globaldata->Groups.size() << endl;
-                       
+                       if ((outputForm == "lt") || (outputForm == "square")) {
+                               //output numSeqs
+                               out << globaldata->Groups.size() << endl;
+                       }
+
                        //make matrix with scores in it
                        vector< vector<float> > dists;  dists.resize(globaldata->Groups.size());
                        for (int i = 0; i < globaldata->Groups.size(); i++) {
@@ -514,7 +545,7 @@ void UnifracWeightedCommand::createPhylipFile() {
                        
                        //flip it so you can print it
                        for (int r=0; r<globaldata->Groups.size(); r++) { 
-                               for (int l = r+1; l < globaldata->Groups.size(); l++) {
+                               for (int l = 0; l < r; l++) {
                                        dists[r][l] = utreeScores[count];
                                        dists[l][r] = utreeScores[count];
                                        count++;
@@ -528,11 +559,30 @@ void UnifracWeightedCommand::createPhylipFile() {
                                if (name.length() < 10) { //pad with spaces to make compatible
                                        while (name.length() < 10) {  name += " ";  }
                                }
-                               out << name << '\t';
                                
-                               //output distances
-                               for (int l = 0; l < r; l++) {   out  << dists[r][l] << '\t';  }
-                               out << endl;
+                               if (outputForm == "lt") {
+                                       out << name << '\t';
+                                       
+                                       //output distances
+                                       for (int l = 0; l < r; l++) {   out  << dists[r][l] << '\t';  }
+                                       out << endl;
+                               }else if (outputForm == "square") {
+                                       out << name << '\t';
+                                       
+                                       //output distances
+                                       for (int l = 0; l < globaldata->Groups.size(); l++) {   out  << dists[r][l] << '\t';  }
+                                       out << endl;
+                               }else{
+                                       //output distances
+                                       for (int l = 0; l < r; l++) {   
+                                               string otherName = globaldata->Groups[l];
+                                               if (otherName.length() < 10) { //pad with spaces to make compatible
+                                                       while (otherName.length() < 10) {  otherName += " ";  }
+                                               }
+                                               
+                                               out  << name << '\t' << otherName << dists[r][l] << endl;  
+                                       }
+                               }
                        }
                        out.close();
                }