]> git.donarmstrong.com Git - mothur.git/blobdiff - removelineagecommand.cpp
removed read.dist, read.otu, read.tree and globaldata. added current to defaults...
[mothur.git] / removelineagecommand.cpp
index 629416a5f01476442b3dd3e8db9cd3fadaca4165..021b56cd02f36aaa88edac90a8773ed1090bc5a0 100644 (file)
 #include "listvector.hpp"
 
 //**********************************************************************************************************************
-vector<string> RemoveLineageCommand::getValidParameters(){     
+vector<string> RemoveLineageCommand::setParameters(){  
        try {
-               string Array[] =  {"fasta","name", "group", "alignreport", "taxon", "dups", "list","taxonomy","outputdir","inputdir"};
-               vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
+               CommandParameter pfasta("fasta", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pfasta);
+               CommandParameter pname("name", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pname);
+               CommandParameter pgroup("group", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(pgroup);
+               CommandParameter plist("list", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(plist);
+               CommandParameter ptaxonomy("taxonomy", "InputTypes", "", "", "none", "FNGLT", "none",false,true); parameters.push_back(ptaxonomy);
+               CommandParameter palignreport("alignreport", "InputTypes", "", "", "none", "FNGLT", "none",false,false); parameters.push_back(palignreport);
+               CommandParameter ptaxon("taxon", "String", "", "", "", "", "",false,true); parameters.push_back(ptaxon);
+               CommandParameter pdups("dups", "Boolean", "", "T", "", "", "",false,false); parameters.push_back(pdups);
+               CommandParameter pinputdir("inputdir", "String", "", "", "", "", "",false,false); parameters.push_back(pinputdir);
+               CommandParameter poutputdir("outputdir", "String", "", "", "", "", "",false,false); parameters.push_back(poutputdir);
+               
+               vector<string> myArray;
+               for (int i = 0; i < parameters.size(); i++) {   myArray.push_back(parameters[i].name);          }
                return myArray;
        }
        catch(exception& e) {
-               m->errorOut(e, "RemoveLineageCommand", "getValidParameters");
+               m->errorOut(e, "RemoveLineageCommand", "setParameters");
                exit(1);
        }
 }
+//**********************************************************************************************************************
+string RemoveLineageCommand::getHelpString(){  
+       try {
+               string helpString = "";
+               helpString += "The remove.lineage command reads a taxonomy file and any of the following file types: fasta, name, group, list or alignreport file.\n";
+               helpString += "It outputs a file containing only the sequences from the taxonomy file that are not from the taxon you requested to be removed.\n";
+               helpString += "The remove.lineage command parameters are taxon, fasta, name, group, list, taxonomy, alignreport and dups.  You must provide taxonomy unless you have a valid current taxonomy file.\n";
+               helpString += "The dups parameter allows you to add the entire line from a name file if you add any name from the line. default=false. \n";
+               helpString += "The taxon parameter allows you to select the taxons you would like to remove, and is required.\n";
+               helpString += "You may enter your taxons with confidence scores, doing so will remove only those sequences that belong to the taxonomy and whose cofidence scores fall below the scores you give.\n";
+               helpString += "If they belong to the taxonomy and have confidences above those you provide the sequence will not be removed.\n";
+               helpString += "The remove.lineage command should be in the following format: remove.lineage(taxonomy=yourTaxonomyFile, taxon=yourTaxons).\n";
+               helpString += "Example remove.lineage(taxonomy=amazon.silva.taxonomy, taxon=Bacteria;Firmicutes;Bacilli;Lactobacillales;).\n";
+               helpString += "Note: If you are running mothur in script mode you must wrap the taxon in ' characters so mothur will ignore the ; in the taxon.\n";
+               helpString += "Example remove.lineage(taxonomy=amazon.silva.taxonomy, taxon='Bacteria;Firmicutes;Bacilli;Lactobacillales;').\n";
+               helpString += "Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n\n";
+               return helpString;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "RemoveLineageCommand", "getHelpString");
+               exit(1);
+       }
+}
+
+
 //**********************************************************************************************************************
 RemoveLineageCommand::RemoveLineageCommand(){  
        try {
-               //initialize outputTypes
+               abort = true; calledHelp = true; 
+               setParameters();
                vector<string> tempOutNames;
                outputTypes["fasta"] = tempOutNames;
                outputTypes["taxonomy"] = tempOutNames;
@@ -41,35 +78,12 @@ RemoveLineageCommand::RemoveLineageCommand(){
        }
 }
 //**********************************************************************************************************************
-vector<string> RemoveLineageCommand::getRequiredParameters(){  
-       try {
-               string Array[] =  {"taxonomy"};
-               vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
-               return myArray;
-       }
-       catch(exception& e) {
-               m->errorOut(e, "RemoveLineageCommand", "getRequiredParameters");
-               exit(1);
-       }
-}
-//**********************************************************************************************************************
-vector<string> RemoveLineageCommand::getRequiredFiles(){       
-       try {
-               vector<string> myArray;
-               return myArray;
-       }
-       catch(exception& e) {
-               m->errorOut(e, "RemoveLineageCommand", "getRequiredFiles");
-               exit(1);
-       }
-}
-//**********************************************************************************************************************
 RemoveLineageCommand::RemoveLineageCommand(string option)  {
        try {
-               abort = false;
+               abort = false; calledHelp = false;   
                                
                //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
@@ -177,10 +191,18 @@ RemoveLineageCommand::RemoveLineageCommand(string option)  {
                        
                        taxfile = validParameter.validFile(parameters, "taxonomy", true);
                        if (taxfile == "not open") { abort = true; }
-                       else if (taxfile == "not found") {  taxfile = ""; m->mothurOut("The taxonomy parameter is required for the get.lineage command."); m->mothurOutEndLine();  abort = true; }
+                       else if (taxfile == "not found") {              
+                               taxfile = m->getTaxonomyFile(); 
+                               if (taxfile != "") { m->mothurOut("Using " + taxfile + " as input file for the taxonomy parameter."); m->mothurOutEndLine(); }
+                               else {  m->mothurOut("You have no current taxonomy file and the taxonomy parameter is required."); m->mothurOutEndLine(); abort = true; }
+                       }
                        
                        string usedDups = "true";
-                       string temp = validParameter.validFile(parameters, "dups", false);      if (temp == "not found") { temp = "false"; usedDups = ""; }
+                       string temp = validParameter.validFile(parameters, "dups", false);      
+                       if (temp == "not found") { 
+                               if (namefile != "") {  temp = "true";                                   }
+                               else                            {  temp = "false"; usedDups = "";       }
+                       }
                        dups = m->isTrue(temp);
                        
                        taxons = validParameter.validFile(parameters, "taxon", false);  
@@ -206,31 +228,10 @@ RemoveLineageCommand::RemoveLineageCommand(string option)  {
 }
 //**********************************************************************************************************************
 
-void RemoveLineageCommand::help(){
-       try {
-               m->mothurOut("The remove.lineage command reads a taxonomy file and any of the following file types: fasta, name, group, list or alignreport file.\n");
-               m->mothurOut("It outputs a file containing only the sequences from the taxonomy file that are not from the taxon you requested to be removed.\n");
-               m->mothurOut("The remove.lineage command parameters are taxon, fasta, name, group, list, taxonomy, alignreport and dups.  You must provide taxonomy and taxon.\n");
-               m->mothurOut("The dups parameter allows you to add the entire line from a name file if you add any name from the line. default=false. \n");
-               m->mothurOut("The taxon parameter allows you to select the taxons you would like to remove.\n");
-               m->mothurOut("The remove.lineage command should be in the following format: remove.lineage(taxonomy=yourTaxonomyFile, taxon=yourTaxons).\n");
-               m->mothurOut("Example remove.lineage(taxonomy=amazon.silva.taxonomy, taxon=Bacteria;Firmicutes;Bacilli;Lactobacillales;).\n");
-               m->mothurOut("Note: If you are running mothur in script mode you must wrap the taxon in ' characters so mothur will ignore the ; in the taxon.\n");
-               m->mothurOut("Example remove.lineage(taxonomy=amazon.silva.taxonomy, taxon='Bacteria;Firmicutes;Bacilli;Lactobacillales;').\n");
-               m->mothurOut("Note: No spaces between parameter labels (i.e. fasta), '=' and parameters (i.e.yourFasta).\n\n");
-       }
-       catch(exception& e) {
-               m->errorOut(e, "RemoveLineageCommand", "help");
-               exit(1);
-       }
-}
-
-//**********************************************************************************************************************
-
 int RemoveLineageCommand::execute(){
        try {
                
-               if (abort == true) { return 0; }
+               if (abort == true) { if (calledHelp) { return 0; }  return 2;   }
                
                if (m->control_pressed) { return 0; }
                
@@ -250,6 +251,33 @@ int RemoveLineageCommand::execute(){
                        m->mothurOut("Output File Names: "); m->mothurOutEndLine();
                        for (int i = 0; i < outputNames.size(); i++) {  m->mothurOut(outputNames[i]); m->mothurOutEndLine();    }
                        m->mothurOutEndLine();
+                       
+                       //set fasta file as new current fastafile
+                       string current = "";
+                       itTypes = outputTypes.find("fasta");
+                       if (itTypes != outputTypes.end()) {
+                               if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setFastaFile(current); }
+                       }
+                       
+                       itTypes = outputTypes.find("name");
+                       if (itTypes != outputTypes.end()) {
+                               if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setNameFile(current); }
+                       }
+                       
+                       itTypes = outputTypes.find("group");
+                       if (itTypes != outputTypes.end()) {
+                               if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setGroupFile(current); }
+                       }
+                       
+                       itTypes = outputTypes.find("list");
+                       if (itTypes != outputTypes.end()) {
+                               if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setListFile(current); }
+                       }
+                       
+                       itTypes = outputTypes.find("taxonomy");
+                       if (itTypes != outputTypes.end()) {
+                               if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setTaxonomyFile(current); }
+                       }
                }
                
                return 0;               
@@ -401,16 +429,8 @@ int RemoveLineageCommand::readName(){
                        in >> secondCol;                        
 
                        vector<string> parsedNames;
-                       //parse second column saving each name
-                       while (secondCol.find_first_of(',') != -1) { 
-                               name = secondCol.substr(0,secondCol.find_first_of(','));
-                               secondCol = secondCol.substr(secondCol.find_first_of(',')+1, secondCol.length());
-                               parsedNames.push_back(name);
-                       }
+                       m->splitAtComma(secondCol, parsedNames);
                        
-                       //get name after last ,
-                       parsedNames.push_back(secondCol);
-
                        vector<string> validSecond;  validSecond.clear();
                        for (int i = 0; i < parsedNames.size(); i++) {
                                if (names.count(parsedNames[i]) == 0) {
@@ -522,6 +542,17 @@ int RemoveLineageCommand::readTax(){
                
                bool wroteSomething = false;
                
+               bool taxonsHasConfidence = false;
+               vector< map<string, float> > searchTaxons;
+               string noConfidenceTaxons = taxons;
+               int hasConPos = taxons.find_first_of('(');
+               if (hasConPos != string::npos) {  
+                       taxonsHasConfidence = true; 
+                       searchTaxons = getTaxons(taxons); 
+                       noConfidenceTaxons = removeConfidences(taxons);
+               }
+               
+               
                while(!in.eof()){
 
                        if (m->control_pressed) { in.close(); out.close(); remove(outputFileName.c_str());  return 0; }
@@ -531,20 +562,99 @@ int RemoveLineageCommand::readTax(){
                        
                        string newtax = tax;
                        
-                       //if the users file contains confidence scores we want to ignore them when searching for the taxons
-                       int hasConfidences = tax.find_first_of('(');
-                       if (hasConfidences != string::npos) { 
-                               newtax = removeConfidences(tax);
+                       //if the users file contains confidence scores we want to ignore them when searching for the taxons, unless the taxon has them
+                       if (!taxonsHasConfidence) {
+                               
+                               int hasConfidences = tax.find_first_of('(');
+                               if (hasConfidences != string::npos) { 
+                                       newtax = removeConfidences(tax);
+                               }
+                               
+                               int pos = newtax.find(taxons);
+                               
+                               if (pos == string::npos) { 
+                                       wroteSomething = true;
+                                       out << name << '\t' << tax << endl;
+                               }else{ //this sequence contains the taxon the user wants to remove
+                                       names.insert(name);
+                               }
+                               
+                       }else{//if taxons has them and you don't them remove taxons
+                               int hasConfidences = tax.find_first_of('(');
+                               if (hasConfidences == string::npos) { 
+                               
+                                       int pos = newtax.find(noConfidenceTaxons);
+                                       
+                                       if (pos == string::npos) { 
+                                               wroteSomething = true;
+                                               out << name << '\t' << tax << endl;
+                                       }else{ //this sequence contains the taxon the user wants to remove
+                                               names.insert(name);
+                                       }
+                               }else { //both have confidences so we want to make sure the users confidences are greater then or equal to the taxons
+                                       //first remove confidences from both and see if the taxonomy exists
+                               
+                                       string noNewTax = tax;
+                                       int hasConfidences = tax.find_first_of('(');
+                                       if (hasConfidences != string::npos) { 
+                                               noNewTax = removeConfidences(tax);
+                                       }
+                                       
+                                       int pos = noNewTax.find(noConfidenceTaxons);
+                                       
+                                       if (pos != string::npos) { //if yes, then are the confidences okay
+                                               
+                                               bool remove = false;
+                                               vector< map<string, float> > usersTaxon = getTaxons(newtax);
+                                               
+                                               //the usersTaxon is most likely longer than the searchTaxons, and searchTaxon[0] may relate to userTaxon[4]
+                                               //we want to "line them up", so we will find the the index where the searchstring starts
+                                               int index = 0;
+                                               for (int i = 0; i < usersTaxon.size(); i++) {
+                                                       
+                                                       if (usersTaxon[i].begin()->first == searchTaxons[0].begin()->first) { 
+                                                               index = i;  
+                                                               int spot = 0;
+                                                               bool goodspot = true;
+                                                               //is this really the start, or are we dealing with a taxon of the same name?
+                                                               while ((spot < searchTaxons.size()) && ((i+spot) < usersTaxon.size())) {
+                                                                       if (usersTaxon[i+spot].begin()->first != searchTaxons[spot].begin()->first) { goodspot = false; break; }
+                                                                       else { spot++; }
+                                                               }
+                                                               
+                                                               if (goodspot) { break; }
+                                                       }
+                                               }
+                                               
+                                               for (int i = 0; i < searchTaxons.size(); i++) {
+                                                       
+                                                       if ((i+index) < usersTaxon.size()) { //just in case, should never be false
+                                                               if (usersTaxon[i+index].begin()->second < searchTaxons[i].begin()->second) { //is the users cutoff less than the search taxons
+                                                                       remove = true;
+                                                                       break;
+                                                               }
+                                                       }else {
+                                                               remove = true;
+                                                               break;
+                                                       }
+                                               }
+                                               
+                                               //passed the test so remove you
+                                               if (remove) {
+                                                       names.insert(name);
+                                               }else {
+                                                       wroteSomething = true;
+                                                       out << name << '\t' << tax << endl;
+                                               }
+                                       }else {
+                                               wroteSomething = true;
+                                               out << name << '\t' << tax << endl;
+                                       }
+                               }
                        }
                        
-                       int pos = newtax.find(taxons);
                        
-                       if (pos == string::npos) { 
-                               wroteSomething = true;
-                               out << name << '\t' << tax << endl;
-                       }else{ //this sequence contains the taxon the user wants to remove
-                               names.insert(name);
-                       }
+                       
                        
                        m->gobble(in);
                }
@@ -563,6 +673,48 @@ int RemoveLineageCommand::readTax(){
        }
 }
 /**************************************************************************************************/
+vector< map<string, float> > RemoveLineageCommand::getTaxons(string tax) {
+       try {
+               
+               vector< map<string, float> > t;
+               string taxon = "";
+               int taxLength = tax.length();
+               for(int i=0;i<taxLength;i++){
+                       if(tax[i] == ';'){
+                               
+                               int openParen = taxon.find_first_of('(');
+                               int closeParen = taxon.find_last_of(')');
+                               
+                               string newtaxon, confidence;
+                               if ((openParen != string::npos) && (closeParen != string::npos)) {
+                                       newtaxon = taxon.substr(0, openParen); //rip off confidence
+                                       confidence = taxon.substr((openParen+1), (closeParen-openParen-1));  
+                               }else{
+                                       newtaxon = taxon;
+                                       confidence = "0";
+                               }
+                               float con = 0;
+                               convert(confidence, con);
+                               
+                               map<string, float> temp;
+                               temp[newtaxon] = con;
+                               t.push_back(temp);
+                               
+                               taxon = "";
+                       }
+                       else{
+                               taxon += tax[i];
+                       }
+               }
+               
+               return t;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "RemoveLineageCommand", "getTaxons");
+               exit(1);
+       }
+}
+/**************************************************************************************************/
 string RemoveLineageCommand::removeConfidences(string tax) {
        try {