/*****************************************************************/
Tree::Tree() {
try {
-
globaldata = GlobalData::getInstance();
- numLeaves = globaldata->gTreemap->getNumSeqs();
+
+ if (globaldata->runParse == true) { parseTreeFile(); globaldata->runParse = false; }
+
+ numLeaves = globaldata->Treenames.size();
numNodes = 2*numLeaves - 1;
tree.resize(numNodes);
for (int i = 0; i < numNodes; i++) {
//initialize leaf nodes
if (i <= (numLeaves-1)) {
- tree[i].setName(globaldata->gTreemap->namesOfSeqs[i]);
- tree[i].setGroup(globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i]));
- //the node knows its index
- tree[i].setIndex(i);
+ tree[i].setName(globaldata->Treenames[i]);
+ tree[i].setGroup(globaldata->gTreemap->getGroup(globaldata->Treenames[i]));
//set pcount and pGroup for groupname to 1.
- tree[i].pcount[globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])] = 1;
- tree[i].pGroups[globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])] = 1;
+ tree[i].pcount[globaldata->gTreemap->getGroup(globaldata->Treenames[i])] = 1;
+ tree[i].pGroups[globaldata->gTreemap->getGroup(globaldata->Treenames[i])] = 1;
//Treemap knows name, group and index to speed up search
- globaldata->gTreemap->setIndex(globaldata->gTreemap->namesOfSeqs[i], i);
+ globaldata->gTreemap->setIndex(globaldata->Treenames[i], i);
//intialize non leaf nodes
}else if (i > (numLeaves-1)) {
tree[i].setName("");
tree[i].setGroup("");
- //the node knows its index
- tree[i].setIndex(i);
}
}
}
}
/*****************************************************************/
-
+Tree::~Tree() {}
+/*****************************************************************/
int Tree::getIndex(string searchName) {
try {
//Treemap knows name, group and index to speed up search
try {
int lc = tree[i].getLChild();
int rc = tree[i].getRChild();
-
+
//set parsimony groups to left child
map<string,int> parsimony = tree[lc].pGroups;
for(it=parsimony.begin();it!=parsimony.end();it++){
if(it->second == 1){
parsimony.erase(it->first);
- it--;
+// it--;
}
}
//set one remaining groups to 1
// p[white] = 1 and p[black] = 1. Now go up a level and merge that with a node who has p[white] = 1
//and you get p[white] = 2, p[black] = 1, but you erase the p[black] because you have a p value higher than 1.
-map<string, int> Tree::mergeUserGroups(int i) {
+map<string, int> Tree::mergeUserGroups(int i, vector<string> g) {
try {
int lc = tree[i].getLChild();
//loop through nodes groups removing the ones the user doesn't want
for (it = tree[lc].pGroups.begin(); it != tree[lc].pGroups.end(); it++) {
- if (inUsersGroups(it->first, globaldata->Groups) != true) { tree[lc].pGroups.erase(it->first); }
+ if (inUsersGroups(it->first, g) != true) { tree[lc].pGroups.erase(it->first); }
}
//loop through nodes groups removing the ones the user doesn't want
for (it = tree[rc].pGroups.begin(); it != tree[rc].pGroups.end(); it++) {
- if (inUsersGroups(it->first, globaldata->Groups) != true) { tree[rc].pGroups.erase(it->first); }
+ if (inUsersGroups(it->first, g) != true) { tree[rc].pGroups.erase(it->first); }
}
//set parsimony groups to left child
maxPars = parsimony[it->first];
}
}
-
+
// this is true if right child had a greater parsimony for a certain group
if(maxPars > 1){
//erase all the groups that are only 1 because you found something with 2.
for(it=parsimony.begin();it!=parsimony.end();it++){
if(it->second == 1){
parsimony.erase(it->first);
- it--;
}
}
- //set one remaining groups to 1
- //so with our above example p[white] = 2 would be left and it would become p[white] = 1
for(it=parsimony.begin();it!=parsimony.end();it++){
parsimony[it->first] = 1;
}
+ }
- }
-
return parsimony;
}
catch(exception& e) {
}
/**************************************************************************************************/
-void Tree::randomLabels() {
+void Tree::randomLabels(vector<string> g) {
try {
- //set up the groups the user wants to include
- setGroups();
-
- for(int i=numLeaves-1;i>=0;i--){
- if(tree[i].pGroups.size() == 0){
- continue;
- }
-
- int escape = 1;
+ for(int i = 0; i < numLeaves; i++){
int z;
-
- while(escape == 1){
- //get random index to switch with
- z = int((float)(i+1) * (float)(rand()) / ((float)RAND_MAX+1.0));
-
- if(tree[z].pGroups.size() != 0){
- escape = 0;
- }
- }
+ //get random index to switch with
+ z = int((float)(i+1) * (float)(rand()) / ((float)RAND_MAX+1.0));
//you only want to randomize the nodes that are from a group the user wants analyzed, so
//if either of the leaf nodes you are about to switch are not in the users groups then you don't want to switch them.
bool treez, treei;
-
- //leaves have only one group so you can just set it to begin()
- it = tree[z].pGroups.begin();
- treez = inUsersGroups(it->first, globaldata->Groups);
-
- it = tree[i].pGroups.begin();
- treei = inUsersGroups(it->first, globaldata->Groups);
+
+ treez = inUsersGroups(tree[z].getGroup(), g);
+ treei = inUsersGroups(tree[i].getGroup(), g);
if ((treez == true) && (treei == true)) {
//switches node i and node z's info.
map<string,int> lib_hold = tree[z].pGroups;
tree[z].pGroups = (tree[i].pGroups);
tree[i].pGroups = (lib_hold);
-
- tree[z].setGroup(tree[z].pGroups.begin()->first);
- tree[i].setGroup(tree[i].pGroups.begin()->first);
-
+
+ string zgroup = tree[z].getGroup();
+ tree[z].setGroup(tree[i].getGroup());
+ tree[i].setGroup(zgroup);
+
+ string zname = tree[z].getName();
+ tree[z].setName(tree[i].getName());
+ tree[i].setName(zname);
+
map<string,int> gcount_hold = tree[z].pcount;
tree[z].pcount = (tree[i].pcount);
tree[i].pcount = (gcount_hold);
}
/**************************************************************************************************/
+void Tree::randomLabels(string groupA, string groupB) {
+ try {
+ int numSeqsA = globaldata->gTreemap->seqsPerGroup[groupA];
+ int numSeqsB = globaldata->gTreemap->seqsPerGroup[groupB];
+
+ vector<string> randomGroups(numSeqsA+numSeqsB, groupA);
+ for(int i=numSeqsA;i<randomGroups.size();i++){
+ randomGroups[i] = groupB;
+ }
+ random_shuffle(randomGroups.begin(), randomGroups.end());
+
+ int randomCounter = 0;
+ for(int i=0;i<numLeaves;i++){
+ if(tree[i].getGroup() == groupA || tree[i].getGroup() == groupB){
+ tree[i].setGroup(randomGroups[randomCounter]);
+ tree[i].pcount.clear();
+ tree[i].pcount[randomGroups[randomCounter]] = 1;
+ tree[i].pGroups.clear();
+ tree[i].pGroups[randomGroups[randomCounter]] = 1;
+ randomCounter++;
+ }
+ }
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function randomLabels. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Tree class function randomLabels. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
+/**************************************************************************************************/
void Tree::randomBlengths() {
try {
for(int i=numNodes-1;i>=0;i--){
}
}
/*************************************************************************************************/
-void Tree::assembleRandomUnifracTree() {
- randomLabels();
+void Tree::assembleRandomUnifracTree(vector<string> g) {
+ randomLabels(g);
+ assembleTree();
+}
+/*************************************************************************************************/
+void Tree::assembleRandomUnifracTree(string groupA, string groupB) {
+ randomLabels(groupA, groupB);
assembleTree();
}
exit(1);
}
}
+/*****************************************************************/
+void Tree::print(ostream& out) {
+ try {
+ int root = findRoot();
+ printBranch(root, out, "branch");
+ out << ";" << endl;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function print. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Tree class function print. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
+/*****************************************************************/
+void Tree::printForBoot(ostream& out) {
+ try {
+ int root = findRoot();
+ printBranch(root, out, "boot");
+ out << ";" << endl;
+ }
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function printForBoot. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Tree class function printForBoot. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
/*****************************************************************/
// This prints out the tree in Newick form.
-void Tree::createNewickFile() {
+void Tree::createNewickFile(string f) {
try {
int root = findRoot();
- filename = getRootName(globaldata->getTreeFile()) + "newick";
+ //filename = getRootName(globaldata->getTreeFile()) + "newick";
+ filename = f;
+
openOutputFile(filename, out);
- printBranch(root);
+ printBranch(root, out, "branch");
// you are at the end of the tree
out << ";" << endl;
+ out.close();
}
catch(exception& e) {
cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function createNewickFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
for (int i = 0; i < numNodes; i++) {
//you found the root
if (tree[i].getParent() == -1) { return i; }
+ //cout << "i = " << i << endl;
+ //cout << "i's parent = " << tree[i].getParent() << endl;
}
return -1;
}
}
/*****************************************************************/
-void Tree::printBranch(int node) {
+void Tree::printBranch(int node, ostream& out, string mode) {
try {
// you are not a leaf
if (tree[node].getLChild() != -1) {
out << "(";
- printBranch(tree[node].getLChild());
+ printBranch(tree[node].getLChild(), out, mode);
out << ",";
- printBranch(tree[node].getRChild());
+ printBranch(tree[node].getRChild(), out, mode);
out << ")";
+ if (mode == "branch") {
+ //if there is a branch length then print it
+ if (tree[node].getBranchLength() != -1) {
+ out << ":" << tree[node].getBranchLength();
+ }
+ }else if (mode == "boot") {
+ //if there is a label then print it
+ if (tree[node].getLabel() != -1) {
+ out << tree[node].getLabel();
+ }
+ }
}else { //you are a leaf
- tree[node].printNode(); //prints out name and branch length
+ out << tree[node].getGroup();
+ if (mode == "branch") {
+ //if there is a branch length then print it
+ if (tree[node].getBranchLength() != -1) {
+ out << ":" << tree[node].getBranchLength();
+ }
+ }else if (mode == "boot") {
+ //if there is a label then print it
+ if (tree[node].getLabel() != -1) {
+ out << tree[node].getLabel();
+ }
+ }
}
}
/*****************************************************************/
-void Tree::setGroups() {
+void Tree::printTree() {
+
+ for(int i=0;i<numNodes;i++){
+ cout << i << '\t';
+ tree[i].printNode();
+ }
+
+}
+
+/*****************************************************************/
+
+void Tree::parseTreeFile() {
+
+ //only takes names from the first tree and assumes that all trees use the same names.
try {
- //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 (globaldata->gTreemap->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);
+ string filename = globaldata->getTreeFile();
+ ifstream filehandle;
+ openInputFile(filename, filehandle);
+ int c, comment;
+ comment = 0;
+
+ //ifyou are not a nexus file
+ if((c = filehandle.peek()) != '#') {
+ while((c = filehandle.peek()) != ';') {
+ while ((c = filehandle.peek()) != ';') {
+ // get past comments
+ if(c == '[') {
+ comment = 1;
+ }
+ if(c == ']'){
+ comment = 0;
+ }
+ if((c == '(') && (comment != 1)){ break; }
+ filehandle.get();
}
+
+ readTreeString(filehandle);
}
-
- //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 < globaldata->gTreemap->namesOfGroups.size(); i++) {
- globaldata->Groups.push_back(globaldata->gTreemap->namesOfGroups[i]);
+ //ifyou are a nexus file
+ }else if((c = filehandle.peek()) == '#') {
+ string holder = "";
+
+ // get past comments
+ while(holder != "translate" && holder != "Translate"){
+ if(holder == "[" || holder == "[!"){
+ comment = 1;
+ }
+ if(holder == "]"){
+ comment = 0;
+ }
+ filehandle >> holder;
+
+ //ifthere is no translate then you must read tree string otherwise use translate to get names
+ if(holder == "tree" && comment != 1){
+ //pass over the "tree rep.6878900 = "
+ while (((c = filehandle.get()) != '(') && ((c = filehandle.peek()) != EOF)) {;}
+
+ if(c == EOF) { break; }
+ filehandle.putback(c); //put back first ( of tree.
+ readTreeString(filehandle);
+ break;
}
}
-
- }else {
- for (int i = 0; i < globaldata->gTreemap->namesOfGroups.size(); i++) {
- globaldata->Groups.push_back(globaldata->gTreemap->namesOfGroups[i]);
+
+ //use nexus translation rather than parsing tree to save time
+ if((holder == "translate") || (holder == "Translate")) {
+
+ string number, name, h;
+ h = ""; // so it enters the loop the first time
+ while((h != ";") && (number != ";")) {
+ filehandle >> number;
+ filehandle >> name;
+
+ //c = , until done with translation then c = ;
+ h = name.substr(name.length()-1, name.length());
+ name.erase(name.end()-1); //erase the comma
+ globaldata->Treenames.push_back(number);
+ }
+ if(number == ";") { globaldata->Treenames.pop_back(); } //in case ';' from translation is on next line instead of next to last name
}
}
+
}
catch(exception& e) {
- cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function setGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function parseTreeFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
exit(1);
}
catch(...) {
- cout << "An unknown error has occurred in the Tree class function setGroups. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ cout << "An unknown error has occurred in the Tree class function parseTreeFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
exit(1);
}
-
}
+/*******************************************************/
-/*****************************************************************/
+/*******************************************************/
+void Tree::readTreeString(ifstream& filehandle) {
+ try {
+ int c;
+ string name; //k
+
+ while((c = filehandle.peek()) != ';') {
+ //ifyou are a name
+ if((c != '(') && (c != ')') && (c != ',') && (c != ':') && (c != '\n') && (c != '\t') && (c != 32)) { //32 is space
+ name = "";
+ c = filehandle.get();
+ // k = c;
+//cout << k << endl;
+ while ((c != '(') && (c != ')') && (c != ',') && (c != ':') && (c != '\n') && (c != 32) && (c != '\t')) {
+ name += c;
+ c = filehandle.get();
+ // k = c;
+//cout << " in name while " << k << endl;
+ }
+
+//cout << "name = " << name << endl;
+ globaldata->Treenames.push_back(name);
+ filehandle.putback(c);
+//k = c;
+//cout << " after putback" << k << endl;
+ }
+
+ if(c == ':') { //read until you reach the end of the branch length
+ while ((c != '(') && (c != ')') && (c != ',') && (c != ';') && (c != '\n') && (c != '\t') && (c != 32)) {
+ c = filehandle.get();
+ // k = c;
+ //cout << " in branch while " << k << endl;
+ }
+ filehandle.putback(c);
+ }
+ c = filehandle.get();
+ if(c == ';') { break; }
+ // k = c;
+//cout << k << endl;
-void Tree::printTree() {
-
- for(int i=0;i<numNodes;i++){
- cout << i << '\t';
- tree[i].printNode();
+ }
}
-
-}
+ catch(exception& e) {
+ cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function parseTreeFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+ catch(...) {
+ cout << "An unknown error has occurred in the Tree class function parseTreeFile. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
+ exit(1);
+ }
+}
-/*****************************************************************/
+/*******************************************************/
+/*******************************************************/
projects / mothur.git / blobdiff