//**********************************************************************************************************************
vector<string> UnifracUnweightedCommand::getValidParameters(){
try {
- string Array[] = {"groups","iters","distance","random", "processors","outputdir","inputdir"};
+ string Array[] = {"groups","iters","distance","random","root", "processors","outputdir","inputdir"};
vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
return myArray;
}
else {
//valid paramters for this command
- string Array[] = {"groups","iters","distance","random", "processors","outputdir","inputdir"};
+ string Array[] = {"groups","iters","distance","random","root", "processors","outputdir","inputdir"};
vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
OptionParser parser(option);
temp = validParameter.validFile(parameters, "random", false); if (temp == "not found") { temp = "f"; }
random = m->isTrue(temp);
+ temp = validParameter.validFile(parameters, "root", false); if (temp == "not found") { temp = "F"; }
+ includeRoot = m->isTrue(temp);
+
temp = validParameter.validFile(parameters, "processors", false); if (temp == "not found"){ temp = "1"; }
convert(temp, processors);
if (numGroups == 1) { numComp++; groupComb.push_back(allGroups); }
- unweighted = new Unweighted(tmap);
+ unweighted = new Unweighted(tmap, includeRoot);
}
void UnifracUnweightedCommand::help(){
try {
m->mothurOut("The unifrac.unweighted command can only be executed after a successful read.tree command.\n");
- m->mothurOut("The unifrac.unweighted command parameters are groups, iters, distance, processors and random. No parameters are required.\n");
+ m->mothurOut("The unifrac.unweighted command parameters are groups, iters, distance, processors, root and random. No parameters are required.\n");
m->mothurOut("The groups parameter allows you to specify which of the groups in your groupfile you would like analyzed. You must enter at least 1 valid group.\n");
m->mothurOut("The group names are separated by dashes. The iters parameter allows you to specify how many random trees you would like compared to your tree.\n");
m->mothurOut("The distance parameter allows you to create a distance file from the results. The default is false. You may set distance to lt, square or column.\n");
m->mothurOut("The random parameter allows you to shut off the comparison to random trees. The default is false, meaning compare don't your trees with randomly generated trees.\n");
+ m->mothurOut("The root parameter allows you to include the entire root in your calculations. The default is false, meaning stop at the root for this comparision instead of the root of the entire tree.\n");
m->mothurOut("The processors parameter allows you to specify the number of processors to use. The default is 1.\n");
m->mothurOut("The unifrac.unweighted command should be in the following format: unifrac.unweighted(groups=yourGroups, iters=yourIters).\n");
m->mothurOut("Example unifrac.unweighted(groups=A-B-C, iters=500).\n");
vector< map<float, float> > rscoreFreq; //map <unweighted score, number of random trees with that score.> -vector entry for each combination.
vector< map<float, float> > rCumul; //map <unweighted score, cumulative percentage of number of random trees with that score or higher.> -vector entry for each combination.
- bool abort, phylip, random;
+ bool abort, phylip, random, includeRoot;
string groups, itersString, outputDir, outputForm;
vector<string> Groups, outputNames; //holds groups to be used
map<string, vector<string> > outputTypes;
//**********************************************************************************************************************
vector<string> UnifracWeightedCommand::getValidParameters(){
try {
- string Array[] = {"groups","iters","distance","random","processors","outputdir","inputdir"};
+ string Array[] = {"groups","iters","distance","random","processors","root","outputdir","inputdir"};
vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
return myArray;
}
else {
//valid paramters for this command
- string Array[] = {"groups","iters","distance","random","processors","outputdir","inputdir"};
+ string Array[] = {"groups","iters","distance","random","processors","root","outputdir","inputdir"};
vector<string> myArray (Array, Array+(sizeof(Array)/sizeof(string)));
OptionParser parser(option);
else { m->mothurOut("Options for distance are: lt, square, or column. Using lt."); m->mothurOutEndLine(); phylip = true; outputForm = "lt"; }
}
- temp = validParameter.validFile(parameters, "random", false); if (temp == "not found") { temp = "F"; }
+ temp = validParameter.validFile(parameters, "random", false); if (temp == "not found") { temp = "F"; }
random = m->isTrue(temp);
+ temp = validParameter.validFile(parameters, "root", false); if (temp == "not found") { temp = "F"; }
+ includeRoot = m->isTrue(temp);
+
temp = validParameter.validFile(parameters, "processors", false); if (temp == "not found"){ temp = "1"; }
convert(temp, processors);
util->setGroups(globaldata->Groups, tmap->namesOfGroups, s, numGroups, "weighted"); //sets the groups the user wants to analyze
util->getCombos(groupComb, globaldata->Groups, numComp);
- weighted = new Weighted(tmap);
+ weighted = new Weighted(tmap, includeRoot);
}
}
void UnifracWeightedCommand::help(){
try {
m->mothurOut("The unifrac.weighted command can only be executed after a successful read.tree command.\n");
- m->mothurOut("The unifrac.weighted command parameters are groups, iters, distance, processors and random. No parameters are required.\n");
+ m->mothurOut("The unifrac.weighted command parameters are groups, iters, distance, processors, root and random. No parameters are required.\n");
m->mothurOut("The groups parameter allows you to specify which of the groups in your groupfile you would like analyzed. You must enter at least 2 valid groups.\n");
m->mothurOut("The group names are separated by dashes. The iters parameter allows you to specify how many random trees you would like compared to your tree.\n");
m->mothurOut("The distance parameter allows you to create a distance file from the results. The default is false.\n");
m->mothurOut("The random parameter allows you to shut off the comparison to random trees. The default is false, meaning don't compare your trees with randomly generated trees.\n");
+ m->mothurOut("The root parameter allows you to include the entire root in your calculations. The default is false, meaning stop at the root for this comparision instead of the root of the entire tree.\n");
m->mothurOut("The processors parameter allows you to specify the number of processors to use. The default is 1.\n");
m->mothurOut("The unifrac.weighted command should be in the following format: unifrac.weighted(groups=yourGroups, iters=yourIters).\n");
m->mothurOut("Example unifrac.weighted(groups=A-B-C, iters=500).\n");
vector< map<float, float> > rCumul; //map <weighted score, cumulative percentage of number of random trees with that score or higher.> -vector entry for each c
map<float, float> validScores; //map contains scores from random
- bool abort, phylip, random;
+ bool abort, phylip, random, includeRoot;
string groups, itersString, outputForm;
vector<string> Groups, outputNames; //holds groups to be used
map<string, vector<string> > outputTypes;
m->mothurOut(namesOfGroupCombos[h][namesOfGroupCombos[h].size()-1]);
m->mothurOut(", skipping."); m->mothurOutEndLine(); results[count] = UW;
}else{
-
+
+ //if including the root this clears rootForGrouping[namesOfGroupCombos[h]]
getRoot(t, nodeBelonging, namesOfGroupCombos[h]);
for(int i=0;i<t->getNumNodes();i++){
map<string, int>::iterator itGroup = t->tree[i].pcount.find(namesOfGroupCombos[h][j]);
if (itGroup != t->tree[i].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
}
-
+
+
//unique calc
if (pcountSize == 0) { }
else if ((t->tree[i].getBranchLength() != -1) && (pcountSize == 1) && (rootForGrouping[namesOfGroupCombos[h]].count(i) == 0)) { //you have a unique branch length and you are not the root
UniqueBL += abs(t->tree[i].getBranchLength());
}
-
+
//total calc
if (pcountSize == 0) { }
else if ((t->tree[i].getBranchLength() != -1) && (pcountSize != 0) && (rootForGrouping[namesOfGroupCombos[h]].count(i) == 0)) { //you have a branch length and you are not the root
totalBL += abs(t->tree[i].getBranchLength());
}
-
}
//cout << UniqueBL << '\t' << totalBL << endl;
UW = (UniqueBL / totalBL);
m->mothurOut(", skipping."); m->mothurOutEndLine(); results[count] = UW;
}else{
+ //if including the root this clears rootForGrouping[namesOfGroupCombos[h]]
getRoot(copyTree, nodeBelonging, namesOfGroupCombos[h]);
for(int i=0;i<copyTree->getNumNodes();i++){
//you are a leaf so get your parent
int index = t->tree[index].getParent();
- //my parent is a potential root
- rootForGrouping[grouping].insert(index);
-
- //while you aren't at root
- while(t->tree[index].getParent() != -1){
-
- if (m->control_pressed) { return 0; }
-
- //am I the root for this grouping? if so I want to stop "early"
- //does my sibling have descendants from the users groups?
- //if so I am not the root
- int parent = t->tree[index].getParent();
- int lc = t->tree[parent].getLChild();
- int rc = t->tree[parent].getRChild();
-
- int sib = lc;
- if (lc == index) { sib = rc; }
+ if (includeRoot) {
+ rootForGrouping[grouping].clear();
+ }else {
+ //my parent is a potential root
+ rootForGrouping[grouping].insert(index);
- map<string, int>::iterator itGroup;
- int pcountSize = 0;
- for (int j = 0; j < grouping.size(); j++) {
- map<string, int>::iterator itGroup = t->tree[sib].pcount.find(grouping[j]);
- if (itGroup != t->tree[sib].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
+ //while you aren't at root
+ while(t->tree[index].getParent() != -1){
+
+ if (m->control_pressed) { return 0; }
+
+ //am I the root for this grouping? if so I want to stop "early"
+ //does my sibling have descendants from the users groups?
+ //if so I am not the root
+ int parent = t->tree[index].getParent();
+ int lc = t->tree[parent].getLChild();
+ int rc = t->tree[parent].getRChild();
+
+ int sib = lc;
+ if (lc == index) { sib = rc; }
+
+ map<string, int>::iterator itGroup;
+ int pcountSize = 0;
+ for (int j = 0; j < grouping.size(); j++) {
+ map<string, int>::iterator itGroup = t->tree[sib].pcount.find(grouping[j]);
+ if (itGroup != t->tree[sib].pcount.end()) { pcountSize++; if (pcountSize > 1) { break; } }
+ }
+
+ //if yes, I am not the root
+ if (pcountSize != 0) {
+ rootForGrouping[grouping].clear();
+ rootForGrouping[grouping].insert(parent);
+ }
+
+ index = parent;
}
- //if yes, I am not the root
- if (pcountSize != 0) {
- rootForGrouping[grouping].clear();
+ //get all nodes above the root to add so we don't add their u values above
+ index = *(rootForGrouping[grouping].begin());
+ while(t->tree[index].getParent() != -1){
+ int parent = t->tree[index].getParent();
rootForGrouping[grouping].insert(parent);
+ //cout << parent << " in root" << endl;
+ index = parent;
}
-
- index = parent;
- }
-
- //get all nodes above the root to add so we don't add their u values above
- index = *(rootForGrouping[grouping].begin());
- while(t->tree[index].getParent() != -1){
- int parent = t->tree[index].getParent();
- rootForGrouping[grouping].insert(parent);
- //cout << parent << " in root" << endl;
- index = parent;
}
return 0;
class Unweighted : public TreeCalculator {
public:
- Unweighted(TreeMap* t) : tmap(t) {};
+ Unweighted(TreeMap* t, bool r) : tmap(t), includeRoot(r) {};
~Unweighted() {};
EstOutput getValues(Tree*, int, string);
EstOutput getValues(Tree*, string, string, int, string);
int processors;
string outputDir;
map< vector<string>, set<int> > rootForGrouping; //maps a grouping combo to the roots for that combo
+ bool includeRoot;
EstOutput driver(Tree*, vector< vector<string> >, int, int);
EstOutput createProcesses(Tree*, vector< vector<string> >);
u -= (double) t->tree[i].pcount[groupB] / (double) tmap->seqsPerGroup[groupB];
}
- //if this is not the root then add it
- if (rootForGrouping[namesOfGroupCombos[h]].count(i) == 0) {
+ if (includeRoot) {
if (t->tree[i].getBranchLength() != -1) {
u = abs(u * t->tree[i].getBranchLength());
WScore[(groupA+groupB)] += u;
}
+ }else {
+ //if this is not the root then add it
+ if (rootForGrouping[namesOfGroupCombos[h]].count(i) == 0) {
+ if (t->tree[i].getBranchLength() != -1) {
+ u = abs(u * t->tree[i].getBranchLength());
+ WScore[(groupA+groupB)] += u;
+ }
+ }
}
}
u -= (double) t->tree[i].pcount[groupB] / (double) tmap->seqsPerGroup[groupB];
}
- //if this is not the root then add it
- if (rootForGrouping[groups].count(i) == 0) {
+ if (includeRoot) {
if (t->tree[i].getBranchLength() != -1) {
u = abs(u * t->tree[i].getBranchLength());
WScore[(groupA+groupB)] += u;
}
+ }else{
+ //if this is not the root then add it
+ if (rootForGrouping[groups].count(i) == 0) {
+ if (t->tree[i].getBranchLength() != -1) {
+ u = abs(u * t->tree[i].getBranchLength());
+ WScore[(groupA+groupB)] += u;
+ }
+ }
}
}
/********************************************************/
while(t->tree[index].getParent() != -1){
if (m->control_pressed) { return sum; }
-
- //am I the root for this grouping? if so I want to stop "early"
- //does my sibling have descendants from the users groups?
+
int parent = t->tree[index].getParent();
- int lc = t->tree[parent].getLChild();
- int rc = t->tree[parent].getRChild();
- int sib = lc;
- if (lc == index) { sib = rc; }
-
- map<string, int>::iterator itGroup;
- int pcountSize = 0;
- itGroup = t->tree[sib].pcount.find(groupA);
- if (itGroup != t->tree[sib].pcount.end()) { pcountSize++; }
- itGroup = t->tree[sib].pcount.find(groupB);
- if (itGroup != t->tree[sib].pcount.end()) { pcountSize++; }
-
- //if yes, I am not the root so add me
- if (pcountSize != 0) {
- if (t->tree[index].getBranchLength() != -1) {
- sum += abs(t->tree[index].getBranchLength()) + tempTotal;
- tempTotal = 0.0;
- }else {
- sum += tempTotal;
- tempTotal = 0.0;
- }
- rootForGrouping[grouping].clear();
- rootForGrouping[grouping].insert(parent);
- }else { //if no, I may be the root so add my br to tempTotal until I am proven innocent
- if (t->tree[index].getBranchLength() != -1) {
- tempTotal += abs(t->tree[index].getBranchLength());
+ if (includeRoot) { //add everyone
+ if(t->tree[index].getBranchLength() != -1){ sum += abs(t->tree[index].getBranchLength()); }
+ }else {
+
+ //am I the root for this grouping? if so I want to stop "early"
+ //does my sibling have descendants from the users groups?
+ int lc = t->tree[parent].getLChild();
+ int rc = t->tree[parent].getRChild();
+
+ int sib = lc;
+ if (lc == index) { sib = rc; }
+
+ map<string, int>::iterator itGroup;
+ int pcountSize = 0;
+ itGroup = t->tree[sib].pcount.find(groupA);
+ if (itGroup != t->tree[sib].pcount.end()) { pcountSize++; }
+ itGroup = t->tree[sib].pcount.find(groupB);
+ if (itGroup != t->tree[sib].pcount.end()) { pcountSize++; }
+
+ //if yes, I am not the root so add me
+ if (pcountSize != 0) {
+ if (t->tree[index].getBranchLength() != -1) {
+ sum += abs(t->tree[index].getBranchLength()) + tempTotal;
+ tempTotal = 0.0;
+ }else {
+ sum += tempTotal;
+ tempTotal = 0.0;
+ }
+ rootForGrouping[grouping].clear();
+ rootForGrouping[grouping].insert(parent);
+ }else { //if no, I may be the root so add my br to tempTotal until I am proven innocent
+ if (t->tree[index].getBranchLength() != -1) {
+ tempTotal += abs(t->tree[index].getBranchLength());
+ }
}
}
class Weighted : public TreeCalculator {
public:
- Weighted(TreeMap* t) : tmap(t) {};
+ Weighted(TreeMap* t, bool r) : tmap(t), includeRoot(r) {};
~Weighted() {};
EstOutput getValues(Tree*, string, string);
int processors;
string outputDir;
map< vector<string>, set<int> > rootForGrouping; //maps a grouping combo to the root for that combo
+ bool includeRoot;
EstOutput driver(Tree*, vector< vector<string> >, int, int);
EstOutput createProcesses(Tree*, vector< vector<string> >);
projects / mothur.git / commitdiff