5 * Created by westcott on 10/13/09.
6 * Copyright 2009 Schloss Lab. All rights reserved.
11 #include "rabundvector.hpp"
12 #include "listvector.hpp"
13 #include "sparsematrix.hpp"
15 /***********************************************************************/
17 HCluster::HCluster(RAbundVector* rav, ListVector* lv) : rabund(rav), list(lv){
20 numSeqs = list->getNumSeqs();
22 //initialize cluster array
23 for (int i = 0; i < numSeqs; i++) {
24 clusterNode temp(1, -1, i);
25 clusterArray.push_back(temp);
30 errorOut(e, "HCluster", "HCluster");
34 /***********************************************************************/
36 void HCluster::clusterBins(){
38 //cout << smallCol << '\t' << smallRow << '\t' << smallDist << '\t' << rabund->get(clusterArray[smallRow].smallChild) << '\t' << rabund->get(clusterArray[smallCol].smallChild);
40 rabund->set(clusterArray[smallCol].smallChild, rabund->get(clusterArray[smallRow].smallChild)+rabund->get(clusterArray[smallCol].smallChild));
41 rabund->set(clusterArray[smallRow].smallChild, 0);
42 rabund->setLabel(toString(smallDist));
44 //cout << '\t' << rabund->get(clusterArray[smallRow].smallChild) << '\t' << rabund->get(clusterArray[smallCol].smallChild) << endl;
47 errorOut(e, "HCluster", "clusterBins");
54 /***********************************************************************/
56 void HCluster::clusterNames(){
58 ///cout << smallCol << '\t' << smallRow << '\t' << smallDist << '\t' << list->get(clusterArray[smallRow].smallChild) << '\t' << list->get(clusterArray[smallCol].smallChild);
60 list->set(clusterArray[smallCol].smallChild, list->get(clusterArray[smallRow].smallChild)+','+list->get(clusterArray[smallCol].smallChild));
61 list->set(clusterArray[smallRow].smallChild, "");
62 list->setLabel(toString(smallDist));
64 //cout << '\t' << list->get(clusterArray[smallRow].smallChild) << '\t' << list->get(clusterArray[smallCol].smallChild) << endl;
68 errorOut(e, "HCluster", "clusterNames");
73 /***********************************************************************/
74 int HCluster::getUpmostParent(int node){
77 while (clusterArray[node].parent != -1) {
78 node = clusterArray[node].parent;
84 errorOut(e, "HCluster", "getUpmostParent");
88 /***********************************************************************/
89 void HCluster::printInfo(){
92 cout << "link table" << endl;
93 for (it = activeLinks.begin(); it!= activeLinks.end(); it++) {
94 cout << it->first << " = " << it->second << endl;
97 for (int i = 0; i < linkTable.size(); i++) {
99 for (it = linkTable[i].begin(); it != linkTable[i].end(); it++) {
100 cout << it->first << '-' << it->second << '\t';
104 cout << endl << "clusterArray" << endl;
106 for (int i = 0; i < clusterArray.size(); i++) {
107 cout << i << '\t' << clusterArray[i].numSeq << '\t' << clusterArray[i].parent << '\t' << clusterArray[i].smallChild << endl;
113 catch(exception& e) {
114 errorOut(e, "HCluster", "getUpmostParent");
118 /***********************************************************************/
119 int HCluster::makeActive() {
123 //cout << "active - here" << endl;
124 it = activeLinks.find(smallRow);
125 it2 = activeLinks.find(smallCol);
127 if ((it == activeLinks.end()) && (it2 == activeLinks.end())) { //both are not active so add them
128 int size = linkTable.size();
129 map<int, int> temp; map<int, int> temp2;
131 //add link to eachother
132 temp[smallRow] = 1; // 1 2
133 temp2[smallCol] = 1; // 1 0 1
135 linkTable.push_back(temp);
136 linkTable.push_back(temp2);
139 activeLinks[smallRow] = size;
140 activeLinks[smallCol] = size+1;
141 //cout << "active - here1" << endl;
142 }else if ((it != activeLinks.end()) && (it2 == activeLinks.end())) { //smallRow is active, smallCol is not
143 int size = linkTable.size();
144 int alreadyActiveRow = it->second;
147 //add link to eachother
148 temp[smallRow] = 1; // 6 2 3 5
149 linkTable.push_back(temp); // 6 0 1 2 0
150 linkTable[alreadyActiveRow][smallCol] = 1; // 2 1 0 1 1
154 activeLinks[smallCol] = size;
155 //cout << "active - here2" << endl;
156 }else if ((it == activeLinks.end()) && (it2 != activeLinks.end())) { //smallCol is active, smallRow is not
157 int size = linkTable.size();
158 int alreadyActiveCol = it2->second;
161 //add link to eachother
162 temp[smallCol] = 1; // 6 2 3 5
163 linkTable.push_back(temp); // 6 0 1 2 0
164 linkTable[alreadyActiveCol][smallRow] = 1; // 2 1 0 1 1
168 activeLinks[smallRow] = size;
169 //cout << "active - here3" << endl;
170 }else { //both are active so add one
171 int row = it->second;
172 int col = it2->second;
173 //cout << row << '\t' << col << endl;
175 linkTable[row][smallCol]++;
176 linkTable[col][smallRow]++;
177 linkValue = linkTable[row][smallCol];
178 //cout << "active - here4" << endl;
183 catch(exception& e) {
184 errorOut(e, "HCluster", "makeActive");
188 /***********************************************************************/
189 void HCluster::updateArrayandLinkTable() {
191 //if cluster was made update clusterArray and linkTable
192 int size = clusterArray.size();
195 clusterNode temp(clusterArray[smallRow].numSeq + clusterArray[smallCol].numSeq, -1, clusterArray[smallCol].smallChild);
196 clusterArray.push_back(temp);
199 clusterArray[smallRow].parent = size;
200 clusterArray[smallCol].parent = size;
202 //update linkTable by merging clustered rows and columns
203 int rowSpot = activeLinks[smallRow];
204 int colSpot = activeLinks[smallCol];
205 //cout << "here" << endl;
207 for (int i = 0; i < linkTable.size(); i++) {
208 //check if they are in map
209 it = linkTable[i].find(smallRow);
210 it2 = linkTable[i].find(smallCol);
212 if ((it!=linkTable[i].end()) && (it2!=linkTable[i].end())) { //they are both there
213 linkTable[i][size] = linkTable[i][smallRow]+linkTable[i][smallCol];
214 linkTable[i].erase(smallCol); //delete col row
215 linkTable[i].erase(smallRow); //delete col row
216 }else if ((it==linkTable[i].end()) && (it2!=linkTable[i].end())) { //only col
217 linkTable[i][size] = linkTable[i][smallCol];
218 linkTable[i].erase(smallCol); //delete col
219 }else if ((it!=linkTable[i].end()) && (it2==linkTable[i].end())) { //only row
220 linkTable[i][size] = linkTable[i][smallRow];
221 linkTable[i].erase(smallRow); //delete col
225 //cout << "here2" << endl;
227 for (it = linkTable[rowSpot].begin(); it != linkTable[rowSpot].end(); it++) {
228 it2 = linkTable[colSpot].find(it->first); //does the col also have this
230 if (it2 == linkTable[colSpot].end()) { //not there so add it
231 linkTable[colSpot][it->first] = it->second;
233 linkTable[colSpot][it->first] = it->second+it2->second;
236 //cout << "here3" << endl;
237 linkTable[colSpot].erase(size);
238 linkTable.erase(linkTable.begin()+rowSpot); //delete row
241 activeLinks.erase(smallRow);
242 activeLinks.erase(smallCol);
243 activeLinks[size] = colSpot;
245 //adjust everybody elses spot since you deleted - time vs. space
246 for (it = activeLinks.begin(); it != activeLinks.end(); it++) {
247 if (it->second > rowSpot) { activeLinks[it->first]--; }
250 //cout << "here4" << endl;
253 catch(exception& e) {
254 errorOut(e, "HCluster", "updateArrayandLinkTable");
258 /***********************************************************************/
259 bool HCluster::update(int row, int col, float distance){
264 smallDist = distance;
265 bool clustered = false;
267 //find upmost parent of row and col
268 smallRow = getUpmostParent(smallRow);
269 smallCol = getUpmostParent(smallCol);
270 //cout << "row = " << row << " smallRow = " << smallRow << " col = " << col << " smallCol = " << smallCol << " dist = " << distance << endl;
272 //are they active in the link table
273 int linkValue = makeActive(); //after this point this nodes info is active in linkTable
275 //cout << "linkValue = " << linkValue << " times = " << (clusterArray[smallRow].numSeq * clusterArray[smallCol].numSeq) << endl;
277 if (linkValue == (clusterArray[smallRow].numSeq * clusterArray[smallCol].numSeq)) {
279 updateArrayandLinkTable();
289 catch(exception& e) {
290 errorOut(e, "HCluster", "update");
296 /***********************************************************************/