2 // getotulabelscommand.cpp
5 // Created by Sarah Westcott on 5/21/12.
6 // Copyright (c) 2012 Schloss Lab. All rights reserved.
9 #include "getotulabelscommand.h"
11 //**********************************************************************************************************************
12 vector<string> GetOtuLabelsCommand::setParameters(){
14 CommandParameter paccnos("accnos", "InputTypes", "", "", "none", "none", "none","",false,true, true); parameters.push_back(paccnos);
15 CommandParameter pconstaxonomy("constaxonomy", "InputTypes", "", "", "none", "FNGLT", "none","constaxonomy",false,false, true); parameters.push_back(pconstaxonomy);
16 CommandParameter plist("list", "InputTypes", "", "", "none", "FNGLT", "none","list",false,false, true); parameters.push_back(plist);
17 CommandParameter pshared("shared", "InputTypes", "", "", "none", "FNGLT", "none","shared",false,false, true); parameters.push_back(pshared);
18 CommandParameter potucorr("otucorr", "InputTypes", "", "", "none", "FNGLT", "none","otucorr",false,false, true); parameters.push_back(potucorr);
19 CommandParameter pcorraxes("corraxes", "InputTypes", "", "", "none", "FNGLT", "none","corraxes",false,false, true); parameters.push_back(pcorraxes);
20 CommandParameter plabel("label", "String", "", "", "", "", "","",false,false); parameters.push_back(plabel);
21 CommandParameter pinputdir("inputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(pinputdir);
22 CommandParameter poutputdir("outputdir", "String", "", "", "", "", "","",false,false); parameters.push_back(poutputdir);
24 vector<string> myArray;
25 for (int i = 0; i < parameters.size(); i++) { myArray.push_back(parameters[i].name); }
29 m->errorOut(e, "GetOtuLabelsCommand", "setParameters");
33 //**********************************************************************************************************************
34 string GetOtuLabelsCommand::getHelpString(){
36 string helpString = "";
37 helpString += "The get.otulabels command can be used to select specific otus with the output from classify.otu, otu.association, or corr.axes commands. It can also be used to select a set of otus from a shared or list file.\n";
38 helpString += "The get.otulabels parameters are: constaxonomy, otucorr, corraxes, shared, list, label and accnos.\n";
39 helpString += "The constaxonomy parameter is used to input the results of the classify.otu command.\n";
40 helpString += "The otucorr parameter is used to input the results of the otu.association command.\n";
41 helpString += "The corraxes parameter is used to input the results of the corr.axes command.\n";
42 helpString += "The label parameter is used to analyze specific labels in your input. \n";
43 helpString += "The get.otulabels commmand should be in the following format: \n";
44 helpString += "get.otulabels(accnos=yourListOfOTULabels, corraxes=yourCorrAxesFile)\n";
48 m->errorOut(e, "GetOtuLabelsCommand", "getHelpString");
52 //**********************************************************************************************************************
53 string GetOtuLabelsCommand::getOutputPattern(string type) {
57 if (type == "constaxonomy") { pattern = "[filename],pick,[extension]"; }
58 else if (type == "otucorr") { pattern = "[filename],pick,[extension]"; }
59 else if (type == "corraxes") { pattern = "[filename],pick,[extension]"; }
60 else if (type == "list") { pattern = "[filename],[distance],pick,[extension]"; }
61 else if (type == "shared") { pattern = "[filename],[distance],pick,[extension]"; }
62 else { m->mothurOut("[ERROR]: No definition for type " + type + " output pattern.\n"); m->control_pressed = true; }
67 m->errorOut(e, "GetOtuLabelsCommand", "getOutputPattern");
72 //**********************************************************************************************************************
73 GetOtuLabelsCommand::GetOtuLabelsCommand(){
75 abort = true; calledHelp = true;
77 vector<string> tempOutNames;
78 outputTypes["constaxonomy"] = tempOutNames;
79 outputTypes["otucorr"] = tempOutNames;
80 outputTypes["corraxes"] = tempOutNames;
81 outputTypes["shared"] = tempOutNames;
82 outputTypes["list"] = tempOutNames;
85 m->errorOut(e, "GetOtuLabelsCommand", "GetOtuLabelsCommand");
89 //**********************************************************************************************************************
90 GetOtuLabelsCommand::GetOtuLabelsCommand(string option) {
92 abort = false; calledHelp = false;
94 //allow user to run help
95 if(option == "help") { help(); abort = true; calledHelp = true; }
96 else if(option == "citation") { citation(); abort = true; calledHelp = true;}
99 //valid paramters for this command
100 vector<string> myArray = setParameters();
102 OptionParser parser(option);
103 map<string,string> parameters = parser.getParameters();
105 ValidParameters validParameter;
106 map<string,string>::iterator it;
107 //check to make sure all parameters are valid for command
108 for (it = parameters.begin(); it != parameters.end(); it++) {
109 if (validParameter.isValidParameter(it->first, myArray, it->second) != true) { abort = true; }
113 //if the user changes the input directory command factory will send this info to us in the output parameter
114 string inputDir = validParameter.validFile(parameters, "inputdir", false);
115 if (inputDir == "not found"){ inputDir = ""; }
118 //edit file types below to include only the types you added as parameters
121 it = parameters.find("constaxonomy");
122 //user has given a template file
123 if(it != parameters.end()){
124 path = m->hasPath(it->second);
125 //if the user has not given a path then, add inputdir. else leave path alone.
126 if (path == "") { parameters["constaxonomy"] = inputDir + it->second; }
129 it = parameters.find("accnos");
130 //user has given a template file
131 if(it != parameters.end()){
132 path = m->hasPath(it->second);
133 //if the user has not given a path then, add inputdir. else leave path alone.
134 if (path == "") { parameters["accnos"] = inputDir + it->second; }
137 it = parameters.find("corraxes");
138 //user has given a template file
139 if(it != parameters.end()){
140 path = m->hasPath(it->second);
141 //if the user has not given a path then, add inputdir. else leave path alone.
142 if (path == "") { parameters["corraxes"] = inputDir + it->second; }
145 it = parameters.find("otucorr");
146 //user has given a template file
147 if(it != parameters.end()){
148 path = m->hasPath(it->second);
149 //if the user has not given a path then, add inputdir. else leave path alone.
150 if (path == "") { parameters["otucorr"] = inputDir + it->second; }
153 it = parameters.find("list");
154 //user has given a template file
155 if(it != parameters.end()){
156 path = m->hasPath(it->second);
157 //if the user has not given a path then, add inputdir. else leave path alone.
158 if (path == "") { parameters["list"] = inputDir + it->second; }
161 it = parameters.find("shared");
162 //user has given a template file
163 if(it != parameters.end()){
164 path = m->hasPath(it->second);
165 //if the user has not given a path then, add inputdir. else leave path alone.
166 if (path == "") { parameters["shared"] = inputDir + it->second; }
170 vector<string> tempOutNames;
171 outputTypes["constaxonomy"] = tempOutNames;
172 outputTypes["otucorr"] = tempOutNames;
173 outputTypes["corraxes"] = tempOutNames;
174 outputTypes["shared"] = tempOutNames;
175 outputTypes["list"] = tempOutNames;
177 //check for parameters
178 accnosfile = validParameter.validFile(parameters, "accnos", true);
179 if (accnosfile == "not open") { abort = true; }
180 else if (accnosfile == "not found") {
181 accnosfile = m->getAccnosFile();
182 if (accnosfile != "") { m->mothurOut("Using " + accnosfile + " as input file for the accnos parameter."); m->mothurOutEndLine(); }
184 m->mothurOut("You have no valid accnos file and accnos is required."); m->mothurOutEndLine();
187 }else { m->setAccnosFile(accnosfile); }
189 constaxonomyfile = validParameter.validFile(parameters, "constaxonomy", true);
190 if (constaxonomyfile == "not open") { constaxonomyfile = ""; abort = true; }
191 else if (constaxonomyfile == "not found") { constaxonomyfile = ""; }
193 corraxesfile = validParameter.validFile(parameters, "corraxes", true);
194 if (corraxesfile == "not open") { corraxesfile = ""; abort = true; }
195 else if (corraxesfile == "not found") { corraxesfile = ""; }
197 otucorrfile = validParameter.validFile(parameters, "otucorr", true);
198 if (otucorrfile == "not open") { otucorrfile = ""; abort = true; }
199 else if (otucorrfile == "not found") { otucorrfile = ""; }
201 listfile = validParameter.validFile(parameters, "list", true);
202 if (listfile == "not open") { listfile = ""; abort = true; }
203 else if (listfile == "not found") { listfile = ""; }
204 else { m->setListFile(listfile); }
206 sharedfile = validParameter.validFile(parameters, "shared", true);
207 if (sharedfile == "not open") { sharedfile = ""; abort = true; }
208 else if (sharedfile == "not found") { sharedfile = ""; }
209 else { m->setSharedFile(sharedfile); }
211 //if the user changes the output directory command factory will send this info to us in the output parameter
212 outputDir = validParameter.validFile(parameters, "outputdir", false); if (outputDir == "not found"){ outputDir = ""; }
214 if ((constaxonomyfile == "") && (corraxesfile == "") && (otucorrfile == "") && (sharedfile == "") && (listfile == "")) { m->mothurOut("You must provide one of the following: constaxonomy, corraxes, otucorr, shared or list."); m->mothurOutEndLine(); abort = true; }
216 if ((sharedfile != "") || (listfile != "")) {
217 label = validParameter.validFile(parameters, "label", false);
218 if (label == "not found") { label = ""; m->mothurOut("You did not provide a label, I will use the first label in your inputfile."); m->mothurOutEndLine(); label=""; }
223 catch(exception& e) {
224 m->errorOut(e, "GetOtuLabelsCommand", "GetOtuLabelsCommand");
228 //**********************************************************************************************************************
230 int GetOtuLabelsCommand::execute(){
233 if (abort == true) { if (calledHelp) { return 0; } return 2; }
235 //get labels you want to keep
236 labels = m->readAccnos(accnosfile);
238 if (m->control_pressed) { return 0; }
240 //read through the correct file and output lines you want to keep
241 if (constaxonomyfile != "") { readClassifyOtu(); }
242 if (corraxesfile != "") { readCorrAxes(); }
243 if (otucorrfile != "") { readOtuAssociation(); }
244 if (listfile != "") { readList(); }
245 if (sharedfile != "") { readShared(); }
247 if (m->control_pressed) { for (int i = 0; i < outputNames.size(); i++) { m->mothurRemove(outputNames[i]); } return 0; }
249 //output files created by command
250 m->mothurOutEndLine();
251 m->mothurOut("Output File Names: "); m->mothurOutEndLine();
252 for (int i = 0; i < outputNames.size(); i++) { m->mothurOut(outputNames[i]); m->mothurOutEndLine(); }
253 m->mothurOutEndLine();
256 itTypes = outputTypes.find("list");
257 if (itTypes != outputTypes.end()) {
258 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setListFile(current); }
261 itTypes = outputTypes.find("shared");
262 if (itTypes != outputTypes.end()) {
263 if ((itTypes->second).size() != 0) { current = (itTypes->second)[0]; m->setSharedFile(current); }
268 catch(exception& e) {
269 m->errorOut(e, "GetOtuLabelsCommand", "execute");
273 //**********************************************************************************************************************
274 int GetOtuLabelsCommand::readClassifyOtu(){
276 string thisOutputDir = outputDir;
277 if (outputDir == "") { thisOutputDir += m->hasPath(constaxonomyfile); }
278 map<string, string> variables;
279 variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(constaxonomyfile));
280 variables["[extension]"] = m->getExtension(constaxonomyfile);
281 string outputFileName = getOutputFileName("constaxonomy", variables);
284 m->openOutputFile(outputFileName, out);
287 m->openInputFile(constaxonomyfile, in);
289 bool wroteSomething = false;
290 int selectedCount = 0;
293 string headers = m->getline(in);
294 out << headers << endl;
298 if (m->control_pressed) { break; }
300 string otu = ""; string tax = "unknown";
303 in >> otu >> size >> tax; m->gobble(in);
305 if (labels.count(otu) != 0) {
306 wroteSomething = true;
309 out << otu << '\t' << size << '\t' << tax << endl;
315 if (wroteSomething == false) { m->mothurOut("Your file does not contain any labels from the .accnos file."); m->mothurOutEndLine(); }
316 outputNames.push_back(outputFileName); outputTypes["constaxonomy"].push_back(outputFileName);
318 m->mothurOut("Selected " + toString(selectedCount) + " otus from your constaxonomy file."); m->mothurOutEndLine();
323 catch(exception& e) {
324 m->errorOut(e, "GetOtuLabelsCommand", "readClassifyOtu");
328 //**********************************************************************************************************************
329 int GetOtuLabelsCommand::readOtuAssociation(){
331 string thisOutputDir = outputDir;
332 if (outputDir == "") { thisOutputDir += m->hasPath(otucorrfile); }
333 map<string, string> variables;
334 variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(otucorrfile));
335 variables["[extension]"] = m->getExtension(otucorrfile);
336 string outputFileName = getOutputFileName("otucorr", variables);
339 m->openOutputFile(outputFileName, out);
342 m->openInputFile(otucorrfile, in);
344 bool wroteSomething = false;
345 int selectedCount = 0;
348 string headers = m->getline(in);
349 out << headers << endl;
353 if (m->control_pressed) { break; }
358 string line = m->getline(in); m->gobble(in);
360 if ((labels.count(otu1) != 0) && (labels.count(otu2) != 0)){
361 wroteSomething = true;
364 out << otu1 << '\t' << otu2 << '\t' << line << endl;
370 if (wroteSomething == false) { m->mothurOut("Your file does not contain any labels from the .accnos file."); m->mothurOutEndLine(); }
371 outputNames.push_back(outputFileName); outputTypes["otucorr"].push_back(outputFileName);
373 m->mothurOut("Selected " + toString(selectedCount) + " lines from your otu.corr file."); m->mothurOutEndLine();
378 catch(exception& e) {
379 m->errorOut(e, "GetOtuLabelsCommand", "readOtuAssociation");
383 //**********************************************************************************************************************
384 int GetOtuLabelsCommand::readCorrAxes(){
386 string thisOutputDir = outputDir;
387 if (outputDir == "") { thisOutputDir += m->hasPath(corraxesfile); }
388 map<string, string> variables;
389 variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(corraxesfile));
390 variables["[extension]"] = m->getExtension(corraxesfile);
391 string outputFileName = getOutputFileName("corraxes", variables);
394 m->openOutputFile(outputFileName, out);
398 m->openInputFile(corraxesfile, in);
400 bool wroteSomething = false;
401 int selectedCount = 0;
404 string headers = m->getline(in);
405 out << headers << endl;
409 if (m->control_pressed) { break; }
413 string line = m->getline(in); m->gobble(in);
415 if (labels.count(otu) != 0) {
416 wroteSomething = true;
419 out << otu << '\t' << line << endl;
425 if (wroteSomething == false) { m->mothurOut("Your file does not contain any labels from the .accnos file."); m->mothurOutEndLine(); }
426 outputNames.push_back(outputFileName); outputTypes["corraxes"].push_back(outputFileName);
428 m->mothurOut("Selected " + toString(selectedCount) + " lines from your corr.axes file."); m->mothurOutEndLine();
433 catch(exception& e) {
434 m->errorOut(e, "GetOtuLabelsCommand", "readCorrAxes");
438 //**********************************************************************************************************************
439 int GetOtuLabelsCommand::readShared(){
444 if (m->control_pressed) { for (int i = 0; i < lookup.size(); i++) { delete lookup[i]; } return 0; }
446 vector<string> newLabels;
448 //create new "filtered" lookup
449 vector<SharedRAbundVector*> newLookup;
450 for (int i = 0; i < lookup.size(); i++) {
451 SharedRAbundVector* temp = new SharedRAbundVector();
452 temp->setLabel(lookup[i]->getLabel());
453 temp->setGroup(lookup[i]->getGroup());
454 newLookup.push_back(temp);
457 bool wroteSomething = false;
459 for (int i = 0; i < lookup[0]->getNumBins(); i++) {
461 if (m->control_pressed) { for (int j = 0; j < newLookup.size(); j++) { delete newLookup[j]; } for (int j = 0; j < lookup.size(); j++) { delete lookup[j]; } return 0; }
463 //is this otu on the list
464 if (labels.count(m->currentBinLabels[i]) != 0) {
465 numSelected++; wroteSomething = true;
466 newLabels.push_back(m->currentBinLabels[i]);
467 for (int j = 0; j < newLookup.size(); j++) { //add this OTU to the new lookup
468 newLookup[j]->push_back(lookup[j]->getAbundance(i), lookup[j]->getGroup());
473 string thisOutputDir = outputDir;
474 if (outputDir == "") { thisOutputDir += m->hasPath(sharedfile); }
475 map<string, string> variables;
476 variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(sharedfile));
477 variables["[extension]"] = m->getExtension(sharedfile);
478 variables["[distance]"] = lookup[0]->getLabel();
479 string outputFileName = getOutputFileName("shared", variables);
481 m->openOutputFile(outputFileName, out);
482 outputTypes["shared"].push_back(outputFileName); outputNames.push_back(outputFileName);
484 for (int j = 0; j < lookup.size(); j++) { delete lookup[j]; }
486 m->currentBinLabels = newLabels;
488 newLookup[0]->printHeaders(out);
490 for (int i = 0; i < newLookup.size(); i++) {
491 out << newLookup[i]->getLabel() << '\t' << newLookup[i]->getGroup() << '\t';
492 newLookup[i]->print(out);
496 for (int j = 0; j < newLookup.size(); j++) { delete newLookup[j]; }
498 if (wroteSomething == false) { m->mothurOut("Your file does not contain any OTUs from the .accnos file."); m->mothurOutEndLine(); }
500 m->mothurOut("Selected " + toString(numSelected) + " OTUs from your shared file."); m->mothurOutEndLine();
504 catch(exception& e) {
505 m->errorOut(e, "GetOtuLabelsCommand", "readShared");
509 //**********************************************************************************************************************
510 int GetOtuLabelsCommand::readList(){
514 if (m->control_pressed) { delete list; return 0;}
517 newList.setLabel(list->getLabel());
518 int selectedCount = 0;
519 bool wroteSomething = false;
520 string snumBins = toString(list->getNumBins());
522 for (int i = 0; i < list->getNumBins(); i++) {
524 if (m->control_pressed) { delete list; return 0;}
526 //create a label for this otu
527 string otuLabel = "Otu";
528 string sbinNumber = toString(i+1);
529 if (sbinNumber.length() < snumBins.length()) {
530 int diff = snumBins.length() - sbinNumber.length();
531 for (int h = 0; h < diff; h++) { otuLabel += "0"; }
533 otuLabel += sbinNumber;
535 if (labels.count(otuLabel) != 0) {
537 newList.push_back(list->get(i));
541 string thisOutputDir = outputDir;
542 if (outputDir == "") { thisOutputDir += m->hasPath(listfile); }
543 map<string, string> variables;
544 variables["[filename]"] = thisOutputDir + m->getRootName(m->getSimpleName(listfile));
545 variables["[extension]"] = m->getExtension(listfile);
546 variables["[distance]"] = list->getLabel();
547 string outputFileName = getOutputFileName("list", variables);
549 m->openOutputFile(outputFileName, out);
552 //print new listvector
553 if (newList.getNumBins() != 0) {
554 wroteSomething = true;
559 if (wroteSomething == false) { m->mothurOut("Your file does not contain any OTUs from the .accnos file."); m->mothurOutEndLine(); }
560 outputNames.push_back(outputFileName); outputTypes["list"].push_back(outputFileName);
562 m->mothurOut("Selected " + toString(selectedCount) + " OTUs from your list file."); m->mothurOutEndLine();
566 catch(exception& e) {
567 m->errorOut(e, "GetOtuLabelsCommand", "readList");
571 //**********************************************************************************************************************
572 int GetOtuLabelsCommand::getListVector(){
574 InputData input(listfile, "list");
575 list = input.getListVector();
576 string lastLabel = list->getLabel();
578 if (label == "") { label = lastLabel; return 0; }
580 //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
581 set<string> labels; labels.insert(label);
582 set<string> processedLabels;
583 set<string> userLabels = labels;
585 //as long as you are not at the end of the file or done wih the lines you want
586 while((list != NULL) && (userLabels.size() != 0)) {
587 if (m->control_pressed) { return 0; }
589 if(labels.count(list->getLabel()) == 1){
590 processedLabels.insert(list->getLabel());
591 userLabels.erase(list->getLabel());
595 if ((m->anyLabelsToProcess(list->getLabel(), userLabels, "") == true) && (processedLabels.count(lastLabel) != 1)) {
596 string saveLabel = list->getLabel();
599 list = input.getListVector(lastLabel);
601 processedLabels.insert(list->getLabel());
602 userLabels.erase(list->getLabel());
604 //restore real lastlabel to save below
605 list->setLabel(saveLabel);
609 lastLabel = list->getLabel();
611 //get next line to process
612 //prevent memory leak
614 list = input.getListVector();
618 if (m->control_pressed) { return 0; }
620 //output error messages about any remaining user labels
621 set<string>::iterator it;
622 bool needToRun = false;
623 for (it = userLabels.begin(); it != userLabels.end(); it++) {
624 m->mothurOut("Your file does not include the label " + *it);
625 if (processedLabels.count(lastLabel) != 1) {
626 m->mothurOut(". I will use " + lastLabel + "."); m->mothurOutEndLine();
629 m->mothurOut(". Please refer to " + lastLabel + "."); m->mothurOutEndLine();
633 //run last label if you need to
634 if (needToRun == true) {
636 list = input.getListVector(lastLabel);
641 catch(exception& e) {
642 m->errorOut(e, "GetOtuLabelsCommand", "getListVector");
646 //**********************************************************************************************************************
647 int GetOtuLabelsCommand::getShared(){
649 InputData input(sharedfile, "sharedfile");
650 lookup = input.getSharedRAbundVectors();
651 string lastLabel = lookup[0]->getLabel();
653 if (label == "") { label = lastLabel; return 0; }
655 //if the users enters label "0.06" and there is no "0.06" in their file use the next lowest label.
656 set<string> labels; labels.insert(label);
657 set<string> processedLabels;
658 set<string> userLabels = labels;
660 //as long as you are not at the end of the file or done wih the lines you want
661 while((lookup[0] != NULL) && (userLabels.size() != 0)) {
662 if (m->control_pressed) { return 0; }
664 if(labels.count(lookup[0]->getLabel()) == 1){
665 processedLabels.insert(lookup[0]->getLabel());
666 userLabels.erase(lookup[0]->getLabel());
670 if ((m->anyLabelsToProcess(lookup[0]->getLabel(), userLabels, "") == true) && (processedLabels.count(lastLabel) != 1)) {
671 string saveLabel = lookup[0]->getLabel();
673 for (int i = 0; i < lookup.size(); i++) { delete lookup[i]; }
674 lookup = input.getSharedRAbundVectors(lastLabel);
676 processedLabels.insert(lookup[0]->getLabel());
677 userLabels.erase(lookup[0]->getLabel());
679 //restore real lastlabel to save below
680 lookup[0]->setLabel(saveLabel);
684 lastLabel = lookup[0]->getLabel();
686 //get next line to process
687 //prevent memory leak
688 for (int i = 0; i < lookup.size(); i++) { delete lookup[i]; }
689 lookup = input.getSharedRAbundVectors();
693 if (m->control_pressed) { return 0; }
695 //output error messages about any remaining user labels
696 set<string>::iterator it;
697 bool needToRun = false;
698 for (it = userLabels.begin(); it != userLabels.end(); it++) {
699 m->mothurOut("Your file does not include the label " + *it);
700 if (processedLabels.count(lastLabel) != 1) {
701 m->mothurOut(". I will use " + lastLabel + "."); m->mothurOutEndLine();
704 m->mothurOut(". Please refer to " + lastLabel + "."); m->mothurOutEndLine();
708 //run last label if you need to
709 if (needToRun == true) {
710 for (int i = 0; i < lookup.size(); i++) { if (lookup[i] != NULL) { delete lookup[i]; } }
711 lookup = input.getSharedRAbundVectors(lastLabel);
716 catch(exception& e) {
717 m->errorOut(e, "GetOtuLabelsCommand", "getShared");
721 //**********************************************************************************************************************