]> git.donarmstrong.com Git - mothur.git/blobdiff - parsimonycommand.cpp
parsimony output changed
[mothur.git] / parsimonycommand.cpp
index e741b5e914770e34b49ac5c244e31b041f8520d2..4529f777358ece0b03a6f53d522a54bd7698ac87 100644 (file)
@@ -22,11 +22,9 @@ ParsimonyCommand::ParsimonyCommand() {
                        T = globaldata->gTree;
                        tmap = globaldata->gTreemap;
                        parsFile = globaldata->getTreeFile() + ".parsimony";
-                       openOutputFile(parsFile, out);
+                       parsFileout = globaldata->getTreeFile() + "temp" + ".parsimony";
                        sumFile = globaldata->getTreeFile() + ".psummary";
                        openOutputFile(sumFile, outSum);
-                       //set users groups to analyze
-                       setGroups();
                }else { //user wants random distribution
                        savetmap = globaldata->gTreemap;
                        getUserInput();
@@ -34,8 +32,11 @@ ParsimonyCommand::ParsimonyCommand() {
                        openOutputFile(parsFile, out);
                }
                
+               //set users groups to analyze
+               setGroups();
                convert(globaldata->getIters(), iters);  //how many random trees to generate
                pars = new Parsimony(tmap);
+               counter = 0;
 
        }
        catch(exception& e) {
@@ -50,20 +51,18 @@ ParsimonyCommand::ParsimonyCommand() {
 /***********************************************************/
 int ParsimonyCommand::execute() {
        try {
-               
+       
+               //get pscore for users tree
+               userData.resize(numComp,0);  //data = AB, AC, BC, ABC.
+               randomData.resize(numComp,0);  //data = AB, AC, BC, ABC.
+               rscoreFreq.resize(numComp);  
+               uscoreFreq.resize(numComp);  
+               rCumul.resize(numComp);  
+               uCumul.resize(numComp);  
+               userTreeScores.resize(numComp);  
+               UScoreSig.resize(numComp); 
                                
                if (randomtree == "") {
-                       //get pscore for users tree
-                       userData.resize(numComp,0);  //data = AB, AC, BC, ABC.
-                       randomData.resize(numComp,0);  //data = AB, AC, BC, ABC.
-                       rscoreFreq.resize(numComp);  
-                       uscoreFreq.resize(numComp);  
-                       rCumul.resize(numComp);  
-                       uCumul.resize(numComp);  
-                       validScores.resize(numComp); 
-                       userTreeScores.resize(numComp);  
-                       UScoreSig.resize(numComp); 
-
                        //get pscores for users trees
                        for (int i = 0; i < T.size(); i++) {
                                cout << "Processing tree " << i+1 << endl;
@@ -71,7 +70,6 @@ int ParsimonyCommand::execute() {
                                
                                //output scores for each combination
                                for(int k = 0; k < numComp; k++) {
-                                       cout << "Tree " << i+1 << " Combination " << groupComb[k] << " parsimony score = " << userData[k] << endl;
                                        //update uscoreFreq
                                        it = uscoreFreq[k].find(userData[k]);
                                        if (it == uscoreFreq[k].end()) {//new score
@@ -79,7 +77,7 @@ int ParsimonyCommand::execute() {
                                        }else{ uscoreFreq[k][userData[k]]++; }
                                        
                                        //add users score to valid scores
-                                       validScores[k][userData[k]] = userData[k];
+                                       validScores[userData[k]] = userData[k];
                                        
                                        //save score for summary file
                                        userTreeScores[k].push_back(userData[k]);
@@ -105,7 +103,7 @@ int ParsimonyCommand::execute() {
                                        }
                        
                                        //add randoms score to validscores
-                                       validScores[r][randomData[r]] = randomData[r];
+                                       validScores[randomData[r]] = randomData[r];
                                }
                                
                                delete randT;
@@ -130,26 +128,25 @@ int ParsimonyCommand::execute() {
                                        }
                        
                                        //add randoms score to validscores
-                                       validScores[r][randomData[r]] = randomData[r];
+                                       validScores[randomData[r]] = randomData[r];
                                }
                                
                                delete randT;
                        }
                }
                
-               float rcumul = 0.0000;
-               float ucumul = 0.0000;
-               
                for(int a = 0; a < numComp; a++) {
-               //this loop fills the cumulative maps and put 0.0000 in the score freq map to make it easier to print.
-                       for (it = validScores[a].begin(); it != validScores[a].end(); it++) { 
+                       float rcumul = 0.0000;
+                       float ucumul = 0.0000;
+                       //this loop fills the cumulative maps and put 0.0000 in the score freq map to make it easier to print.
+                       for (it = validScores.begin(); it != validScores.end(); it++) { 
                                if (randomtree == "") {
                                        it2 = uscoreFreq[a].find(it->first);
                                        //user data has that score 
                                        if (it2 != uscoreFreq[a].end()) { uscoreFreq[a][it->first] /= T.size(); ucumul+= it2->second;  }
                                        else { uscoreFreq[a][it->first] = 0.0000; } //no user trees with that score
                                        //make uCumul map
-                                       uCumul[a][it->first] = ucumul-a;
+                                       uCumul[a][it->first] = ucumul;
                                }
                        
                                //make rscoreFreq map and rCumul
@@ -157,7 +154,7 @@ int ParsimonyCommand::execute() {
                                //get percentage of random trees with that info
                                if (it2 != rscoreFreq[a].end()) {  rscoreFreq[a][it->first] /= iters; rcumul+= it2->second;  }
                                else { rscoreFreq[a][it->first] = 0.0000; } //no random trees with that score
-                               rCumul[a][it->first] = rcumul-a;
+                               rCumul[a][it->first] = rcumul;
                        }
                        
                        //find the signifigance of each user trees score when compared to the random trees and save for printing the summary file
@@ -167,7 +164,7 @@ int ParsimonyCommand::execute() {
                }
                
                printParsimonyFile();
-               printUSummaryFile();
+               if (randomtree == "") { printUSummaryFile(); }
                
                //reset globaldata's treemap if you just did random distrib
                if (randomtree != "") { globaldata->gTreemap = savetmap; }
@@ -193,28 +190,29 @@ int ParsimonyCommand::execute() {
 /***********************************************************/
 void ParsimonyCommand::printParsimonyFile() {
        try {
-               //column headers
-               if (randomtree == "") {
-                       out << "Comb" << '\t' << "Score" << '\t' << "UserFreq" << '\t' << "UserCumul" << '\t' << "RandFreq" << '\t' << "RandCumul" << endl;
-               }else {
-                       out << "Comb" << '\t' << "Score" << '\t' << "RandFreq" << '\t' << "RandCumul" << endl;
-               }
+               vector<double> data;
                
                //format output
                out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
-               
+
                for(int a = 0; a < numComp; a++) {
+                       initFile(groupComb[a]);
                        //print each line
-                       for (it = validScores[a].begin(); it != validScores[a].end(); it++) { 
+                       for (it = validScores.begin(); it != validScores.end(); it++) { 
                                if (randomtree == "") {
-                                       out << setprecision(6)  << groupComb[a] << '\t' << it->first << '\t' << '\t'<< uscoreFreq[a][it->first] << '\t' << uCumul[a][it->first] << '\t' << rscoreFreq[a][it->first] << '\t' << rCumul[a][it->first] << endl
+                                       data.push_back(it->first);  data.push_back(uscoreFreq[a][it->first]); data.push_back(uCumul[a][it->first]); data.push_back(rscoreFreq[a][it->first]); data.push_back(rCumul[a][it->first])
                                }else{
-                                       out << setprecision(6) << groupComb[a] << '\t' << it->first << '\t' << '\t' << rscoreFreq[a][it->first] << '\t' << rCumul[a][it->first] << endl;        
+                                       data.push_back(it->first);  data.push_back(rscoreFreq[a][it->first]); data.push_back(rCumul[a][it->first]); 
                                }
+                               output(data);
+                               data.clear();
                        } 
+                       resetFile();
                }
-               out.close();
                
+               out.close();
+               inFile.close();
+               remove(parsFileout.c_str());
        }
        catch(exception& e) {
                cout << "Standard Error: " << e.what() << " has occurred in the ParsimonyCommand class Function printParsimonyFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
@@ -239,6 +237,7 @@ void ParsimonyCommand::printUSummaryFile() {
                for (int i = 0; i< T.size(); i++) {
                        for(int a = 0; a < numComp; a++) {
                                outSum << setprecision(6) << i+1 << '\t' << groupComb[a] << '\t' << '\t' << userTreeScores[a][i] << '\t' << UScoreSig[a][i] << endl;
+                               cout << setprecision(6) << i+1 << '\t' << groupComb[a] << '\t' << '\t' << userTreeScores[a][i] << '\t' << UScoreSig[a][i] << endl;
                        }
                }
                
@@ -274,6 +273,7 @@ void ParsimonyCommand::getUserInput() {
                                
                        //set tmaps seqsPerGroup
                        tmap->seqsPerGroup[toString(i)] = num;
+                       tmap->namesOfGroups.push_back(toString(i));
                        
                        //set tmaps namesOfSeqs
                        for (int j = 0; j < num; j++) {
@@ -357,8 +357,10 @@ void ParsimonyCommand::setGroups() {
                }
                
                //ABC
-               groupComb.push_back(allGroups);
-               numComp++;
+               if (numComp != 1) {
+                       groupComb.push_back(allGroups);
+                       numComp++;
+               }
                
        }
        catch(exception& e) {
@@ -373,4 +375,100 @@ void ParsimonyCommand::setGroups() {
 }
 /*****************************************************************/
 
+void ParsimonyCommand::initFile(string label){
+       try {
+               if(counter != 0){
+                       openOutputFile(parsFileout, out);
+                       openInputFile(parsFile, inFile);
+
+                       string inputBuffer;
+                       getline(inFile, inputBuffer);
+                       
+                       if (randomtree == "") {
+                               out <<  inputBuffer << '\t' << label + "Score" << '\t' << label + "UserFreq" << '\t' << label + "UserCumul" << '\t' << label + "RandFreq" << '\t' << label + "RandCumul" << endl;
+                       }else {
+                               out <<  inputBuffer << '\t' << label + "Score" << '\t' << label + "RandFreq" << '\t' << label + "RandCumul" << endl;
+                       }
+               }else{
+                       openOutputFile(parsFileout, out);
+                       //column headers
+                       if (randomtree == "") {
+                               out << label + "Score" << '\t' << label + "UserFreq" << '\t' << label + "UserCumul" << '\t' << label + "RandFreq" << '\t' << label + "RandCumul" << endl;
+                       }else {
+                               out << label + "Score" << '\t' << label + "RandFreq" << '\t' << label + "RandCumul" << endl;
+                       }
+               }
+
+               out.setf(ios::fixed, ios::floatfield);
+               out.setf(ios::showpoint);
+       }
+       catch(exception& e) {
+               cout << "Standard Error: " << e.what() << " has occurred in the ParsimonyCommand class Function initFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }
+       catch(...) {
+               cout << "An unknown error has occurred in the ParsimonyCommand class function initFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }
+}
+
+/***********************************************************************/
+
+void ParsimonyCommand::output(vector<double> data){
+       try {
+               if(counter != 0){               
+                       string inputBuffer;
+                       getline(inFile, inputBuffer);
+               
+                       if (randomtree == "") {
+                               out << inputBuffer << '\t' << setprecision(6) << data[0] << '\t' << data[1] << '\t' << data[2] << '\t' << data[3] << '\t' << data[4] << endl;
+                       }else{
+                               out << inputBuffer << '\t' << setprecision(6) << data[0] << '\t' << data[1] << '\t' << data[2] << endl;
+                       }
+               }
+               else{
+                       if (randomtree == "") {
+                               out << setprecision(6) << data[0] << '\t' << data[1] << '\t' << data[2] << '\t' << data[3] << '\t' << data[4] << endl;
+                       }else{
+                               out << setprecision(6) << data[0] << '\t' << data[1] << '\t' << data[2] << endl;
+                       }
+               }
+
+       }
+       catch(exception& e) {
+               cout << "Standard Error: " << e.what() << " has occurred in the ParsimonyCommand class Function output. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }
+       catch(...) {
+               cout << "An unknown error has occurred in the ParsimonyCommand class function output. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }
+}
+
+/***********************************************************************/
+
+void ParsimonyCommand::resetFile(){
+       try {
+               if(counter != 0){
+                       out.close();
+                       inFile.close();
+               }
+               else{
+                       out.close();
+               }
+               counter = 1;
+               
+               remove(parsFile.c_str());
+               rename(parsFileout.c_str(), parsFile.c_str());
+       }
+       catch(exception& e) {
+               cout << "Standard Error: " << e.what() << " has occurred in the ParsimonyCommand class Function resetFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }
+       catch(...) {
+               cout << "An unknown error has occurred in the ParsimonyCommand class function resetFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+               exit(1);
+       }       
+}
+