X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=unifracunweightedcommand.cpp;h=bf2129c8b2bbc2e8e727c88bea45049464947bad;hb=692e0c1f69a78b568dc85cbdcea9fb6c189e2e6c;hp=8adff2d896f811d46e80b806acf45261d759c69b;hpb=d037597badc8d18e235c59f0c1114180edb7f98f;p=mothur.git diff --git a/unifracunweightedcommand.cpp b/unifracunweightedcommand.cpp index 8adff2d..bf2129c 100644 --- a/unifracunweightedcommand.cpp +++ b/unifracunweightedcommand.cpp @@ -23,6 +23,30 @@ UnifracUnweightedCommand::UnifracUnweightedCommand() { distFile = globaldata->getTreeFile() + ".uwdistrib"; openOutputFile(distFile, outDist); + //if the user has not entered specific groups to analyze then do them all + if (globaldata->Groups.size() != 0) { + //check that groups are valid + for (int i = 0; i < globaldata->Groups.size(); i++) { + if (tmap->isValidGroup(globaldata->Groups[i]) != true) { + cout << globaldata->Groups[i] << " is not a valid group, and will be disregarded." << endl; + // erase the invalid group from globaldata->Groups + globaldata->Groups.erase (globaldata->Groups.begin()+i); + } + } + + //if the user only entered invalid groups + if (globaldata->Groups.size() == 0) { + cout << "When using the groups parameter you must have at least 1 valid group. I will run the command using all the groups in your groupfile." << endl; + for (int i = 0; i < tmap->namesOfGroups.size(); i++) { + globaldata->Groups.push_back(tmap->namesOfGroups[i]); + } + } + }else { + for (int i = 0; i < tmap->namesOfGroups.size(); i++) { + globaldata->Groups.push_back(tmap->namesOfGroups[i]); + } + } + convert(globaldata->getIters(), iters); //how many random trees to generate unweighted = new Unweighted(tmap); @@ -93,49 +117,52 @@ int UnifracUnweightedCommand::execute() { outDist << i+1 << '\t' << '\t'<< j+1 << '\t' << '\t' << randomData[0] << endl; } + saveRandomScores(); //save all random scores for unweighted file + //find the signifigance of the score - float rcumul = 0.0000; + float rcumul = 1.0000; for (it = rscoreFreq.begin(); it != rscoreFreq.end(); it++) { + rCumul[it->first] = rcumul; //get percentage of random trees with that info rscoreFreq[it->first] /= iters; - rcumul+= it->second; - rCumul[it->first] = rcumul; + rcumul-= it->second; + } //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; + float ucumul = 1.0000; + 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++) { 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; - + //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 rscoreFreq map and rCumul it2 = totalrscoreFreq.find(it->first); + rCumul[it->first] = rcumul; //get percentage of random trees with that info - if (it2 != totalrscoreFreq.end()) { totalrscoreFreq[it->first] /= (iters*T.size()); rcumul+= it2->second; } + 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; + } printUnweightedFile(); printUWSummaryFile(); - //reset randomTree parameter to 0 - globaldata->setRandomTree("0"); + //reset groups parameter + globaldata->Groups.clear(); delete randT; @@ -207,15 +234,14 @@ void UnifracUnweightedCommand::printUWSummaryFile() { /***********************************************************/ 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; + totalrscoreFreq[it->first] += rscoreFreq[it->first]; }else{ //its a new score - totalrscoreFreq[it->first] = 1; + totalrscoreFreq[it->first] = rscoreFreq[it->first]; } } }