]> git.donarmstrong.com Git - mothur.git/blobdiff - unifracunweightedcommand.cpp
fixed bug in heatmap
[mothur.git] / unifracunweightedcommand.cpp
index 8adff2d896f811d46e80b806acf45261d759c69b..fcff253fc4812bafea60669dcb07c301dfa8de62 100644 (file)
@@ -16,13 +16,16 @@ UnifracUnweightedCommand::UnifracUnweightedCommand() {
                
                T = globaldata->gTree;
                tmap = globaldata->gTreemap;
-               unweightedFile = globaldata->getTreeFile() + ".unweighted";
-               openOutputFile(unweightedFile, out);
                sumFile = globaldata->getTreeFile() + ".uwsummary";
                openOutputFile(sumFile, outSum);
-               distFile = globaldata->getTreeFile() + ".uwdistrib";
-               openOutputFile(distFile, outDist);
 
+               util = new SharedUtil();
+               util->setGroups(globaldata->Groups, tmap->namesOfGroups, allGroups, numGroups, "unweighted");   //sets the groups the user wants to analyze
+               util->getCombos(groupComb, globaldata->Groups, numComp);
+               globaldata->setGroups("");
+               
+               if (numGroups == 1) { numComp++; groupComb.push_back(allGroups); }
+                               
                convert(globaldata->getIters(), iters);  //how many random trees to generate
                unweighted = new Unweighted(tmap);
 
@@ -39,105 +42,85 @@ UnifracUnweightedCommand::UnifracUnweightedCommand() {
 /***********************************************************/
 int UnifracUnweightedCommand::execute() {
        try {
+
+               userData.resize(numComp,0);  //data[0] = unweightedscore 
+               randomData.resize(numComp,0); //data[0] = unweightedscore
+               //create new tree with same num nodes and leaves as users
                
-               //get unweighted for users tree
-               userData.resize(1,0);  //data[0] = unweightedscore 
-               randomData.resize(1,0); //data[0] = unweightedscore
-               
-               //format output
-               outDist.setf(ios::fixed, ios::floatfield); outDist.setf(ios::showpoint);
-               outDist << "Tree#" << '\t' << "Iter" << '\t' << "UWScore" << endl;
+               outSum << "Tree#" << '\t' << "Groups" << '\t'  <<  "UWScore" <<'\t' << "UWSig" <<  endl;
+               cout << "Tree#" << '\t' << "Groups" << '\t'  <<  "UWScore" << '\t' << "UWSig" <<  endl;
                
-               //create new tree with same num nodes and leaves as users
-               randT = new Tree();
-                       
                //get pscores for users trees
                for (int i = 0; i < T.size(); i++) {
-                       cout << "Processing tree " << i+1 << endl;
-                       userData = unweighted->getValues(T[i]);  //userData[0] = unweightedscore
-                       
-                       //update uscoreFreq
-                       it = uscoreFreq.find(userData[0]);
-                       if (it == uscoreFreq.end()) {//new score
-                               uscoreFreq[userData[0]] = 1;
-                       }else{ uscoreFreq[userData[0]]++; }
+                       counter = 0;
                        
-                       //add users score to valid scores
-                       validScores[userData[0]] = userData[0];
+                       output = new ColumnFile(globaldata->getTreeFile()  + toString(i+1) + ".unweighted");
                        
-                       //saves users score
-                       utreeScores.push_back(userData[0]);
+                       //get unweighted for users tree
+                       rscoreFreq.resize(numComp);  
+                       rCumul.resize(numComp);  
+                       utreeScores.resize(numComp);  
+                       UWScoreSig.resize(numComp); 
+
+                       userData = unweighted->getValues(T[i]);  //userData[0] = unweightedscore
                        
-                       //copy T[i]'s info.
-                       randT->getCopy(T[i]); 
+                       //output scores for each combination
+                       for(int k = 0; k < numComp; k++) {
+                               //saves users score
+                               utreeScores[k].push_back(userData[k]);
+
+                       }
                        
                        //get unweighted scores for random trees
                        for (int j = 0; j < iters; j++) {
-                               //create a random tree with same topology as T[i], but different labels
-                               randT->assembleRandomUnifracTree();
-                               //get pscore of random tree
-                               randomData = unweighted->getValues(randT);
-                       
-                               //add trees unweighted score to map of scores
-                               it2 = rscoreFreq.find(randomData[0]);
-                               if (it2 != rscoreFreq.end()) {//already have that score
-                                       rscoreFreq[randomData[0]]++;
-                               }else{//first time we have seen this score
-                                       rscoreFreq[randomData[0]] = 1;
-                               }
+                               //we need a different getValues because when we swap the labels we only want to swap those in each parwise comparison
+                               randomData = unweighted->getValues(T[i], "", "");
+                               
+                               for(int k = 0; k < numComp; k++) {      
+                                       //add trees unweighted score to map of scores
+                                       it2 = rscoreFreq[k].find(randomData[k]);
+                                       if (it2 != rscoreFreq[k].end()) {//already have that score
+                                               rscoreFreq[k][randomData[k]]++;
+                                       }else{//first time we have seen this score
+                                               rscoreFreq[k][randomData[k]] = 1;
+                                       }
                                
-                               //add randoms score to validscores
-                               validScores[randomData[0]] = randomData[0];
+                                       //add randoms score to validscores
+                                       validScores[randomData[k]] = randomData[k];
+                               }
                                
-                               //output info to uwdistrib file
-                               outDist << i+1 << '\t' << '\t'<< j+1 << '\t' << '\t' << randomData[0] << endl;
                        }
-                       
-                       //find the signifigance of the score
-                       float rcumul = 0.0000;
-                       for (it = rscoreFreq.begin(); it != rscoreFreq.end(); it++) { 
-                               //get percentage of random trees with that info
-                               rscoreFreq[it->first] /= iters; 
-                               rcumul+= it->second;  
-                               rCumul[it->first] = rcumul;
+               
+                       for(int a = 0; a < numComp; a++) {
+                               float rcumul = 1.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++) { 
+                                       //make rscoreFreq map and rCumul
+                                       it2 = rscoreFreq[a].find(it->first);
+                                       rCumul[a][it->first] = rcumul;
+                                       //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
+                               }
+                               UWScoreSig[a].push_back(rCumul[a][userData[a]]);
                        }
-                       
-                       //save the signifigance of the users score for printing later
-                       UWScoreSig.push_back(rCumul[userData[0]]);
-                       
-                       saveRandomScores(); //save all random scores for unweighted file
                
-                       //clear random data
-                       rscoreFreq.clear();  //you clear this because in the summary file you want the unweighted signifinance to be relative to these 1000 trees.
-                       rCumul.clear();
-               }
                
-               float ucumul = 0.0000;
-               float rcumul = 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++) { 
-                       it2 = uscoreFreq.find(it->first);
-                       //user data has that score 
-                       if (it2 != uscoreFreq.end()) { uscoreFreq[it->first] /= T.size(); ucumul+= it2->second;  }
-                       else { uscoreFreq[it->first] = 0.0000; } //no user trees with that score
-                       //make uCumul map
-                       uCumul[it->first] = ucumul;
+               
+                       printUnweightedFile();
+                       printUWSummaryFile(i);
                        
-                       //make rscoreFreq map and rCumul
-                       it2 = totalrscoreFreq.find(it->first);
-                       //get percentage of random trees with that info
-                       if (it2 != totalrscoreFreq.end()) {  totalrscoreFreq[it->first] /= (iters*T.size()); rcumul+= it2->second;  }
-                       else { totalrscoreFreq[it->first] = 0.0000; } //no random trees with that score
-                       rCumul[it->first] = rcumul;
+                       delete output;
+                       rscoreFreq.clear(); 
+                       rCumul.clear();  
+                       validScores.clear(); 
+                       utreeScores.clear();  
+                       UWScoreSig.clear(); 
                }
                
-               printUnweightedFile();
-               printUWSummaryFile();
-               
-               //reset randomTree parameter to 0
-               globaldata->setRandomTree("0");
-               
-               delete randT;
+               //reset groups parameter
+               globaldata->Groups.clear(); 
+               outSum.close();
                
                return 0;
                
@@ -154,20 +137,20 @@ int UnifracUnweightedCommand::execute() {
 /***********************************************************/
 void UnifracUnweightedCommand::printUnweightedFile() {
        try {
-               //column headers
-               
-               out << "Score" << '\t' << "UserFreq" << '\t' << "UserCumul" << '\t' << "RandFreq" << '\t' << "RandCumul" << endl;
-                               
-               //format output
-               out.setf(ios::fixed, ios::floatfield); out.setf(ios::showpoint);
-               
-               //print each line
-               for (it = validScores.begin(); it != validScores.end(); it++) { 
-                       out << setprecision(6) << it->first << '\t' << '\t' << uscoreFreq[it->first] << '\t' << uCumul[it->first] << '\t' << totalrscoreFreq[it->first] << '\t' << rCumul[it->first] << endl; 
-               } 
-               
-               out.close();
-               
+               vector<double> data;
+               vector<string> tags;
+               tags.push_back("Score"); tags.push_back("RandFreq"); tags.push_back("RandCumul");
+               
+               for(int a = 0; a < numComp; a++) {
+                       output->initFile(groupComb[a], tags);
+                       //print each line
+                       for (it = validScores.begin(); it != validScores.end(); it++) { 
+                               data.push_back(it->first);  data.push_back(rscoreFreq[a][it->first]); data.push_back(rCumul[a][it->first]); 
+                               output->output(data);
+                               data.clear();
+                       } 
+                       output->resetFile();
+               }
        }
        catch(exception& e) {
                cout << "Standard Error: " << e.what() << " has occurred in the UnifracUnweightedCommand class Function printUnweightedFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
@@ -180,20 +163,27 @@ void UnifracUnweightedCommand::printUnweightedFile() {
 }
 
 /***********************************************************/
-void UnifracUnweightedCommand::printUWSummaryFile() {
+void UnifracUnweightedCommand::printUWSummaryFile(int i) {
        try {
-               //column headers
-               outSum << "Tree#" << '\t'  <<  "UWScore" << '\t' << '\t' << "UWSig" <<  endl;
-               
+                               
                //format output
                outSum.setf(ios::fixed, ios::floatfield); outSum.setf(ios::showpoint);
-               
+                       
                //print each line
-               for (int i = 0; i< T.size(); i++) {
-                       outSum << setprecision(6) << i+1 << '\t' << '\t' << utreeScores[i] << '\t' << UWScoreSig[i] << endl; 
+
+               for(int a = 0; a < numComp; a++) {
+                       outSum << i+1 << '\t';
+                       cout << i+1 << '\t';
+                       
+                       if (UWScoreSig[a][0] > (1/(float)iters)) {
+                               outSum << setprecision(6) << groupComb[a]  << '\t' << utreeScores[a][0] << '\t' << setprecision(globaldata->getIters().length()) << UWScoreSig[a][0] << endl;
+                               cout << setprecision(6)  << groupComb[a]  << '\t' << utreeScores[a][0] << '\t' << setprecision(globaldata->getIters().length()) << UWScoreSig[a][0] << endl; 
+                       }else {
+                               outSum << setprecision(6) << groupComb[a]  << '\t' << utreeScores[a][0] << '\t' << setprecision(globaldata->getIters().length()) << "<" << (1/float(iters)) << endl;
+                               cout << setprecision(6)  << groupComb[a]  << '\t' << utreeScores[a][0] << '\t' << setprecision(globaldata->getIters().length()) << "<" << (1/float(iters)) << endl; 
+                       }
                }
                
-               outSum.close();
        }
        catch(exception& e) {
                cout << "Standard Error: " << e.what() << " has occurred in the UnifracUnweightedCommand class Function printUWSummaryFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
@@ -204,29 +194,7 @@ void UnifracUnweightedCommand::printUWSummaryFile() {
                exit(1);
        }
 }
+
 /***********************************************************/
-void UnifracUnweightedCommand::saveRandomScores() {
-       try {
-               //update total map with new random scores
-               for (it = rscoreFreq.begin(); it != rscoreFreq.end(); it++) { 
-                       //does this score already exist in the total map
-                       it2 = totalrscoreFreq.find(it->first);
-                       //if yes then add them
-                       if (it2 != totalrscoreFreq.end()) { 
-                               it2->second += it->second;
-                       }else{ //its a new score
-                               totalrscoreFreq[it->first] = 1;
-                       }
-               }
-       }
-       catch(exception& e) {
-               cout << "Standard Error: " << e.what() << " has occurred in the UnifracUnweightedCommand class Function saveRandomScores. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
-               exit(1);
-       }
-       catch(...) {
-               cout << "An unknown error has occurred in the UnifracUnweightedCommand class function saveRandomScores. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
-               exit(1);
-       }
-}
 
-/***********************************************************/
\ No newline at end of file
+