heatmap.cpp

   1 /*
   2  *  heatmap.cpp
   3  *  Mothur
   4  *
   5  *  Created by Sarah Westcott on 3/25/09.
   6  *  Copyright 2009 Schloss Lab UMASS Amherst. All rights reserved.
   7  *
   8  */
   9
  10 #include "heatmap.h"
  11
  12 //**********************************************************************************************************************
  13 HeatMap::HeatMap(string sort, string scale, int num, int fsize, string dir, string i){
  14         try {
  15                 m = MothurOut::getInstance();
  16 //              format = globaldata->getFormat();
  17                 sorted = sort;
  18                 scaler = scale;
  19                 outputDir = dir;
  20                 numOTU = num;
  21                 fontSize = fsize;
  22                 inputfile = i;
  23         }
  24         catch(exception& e) {
  25                 m->errorOut(e, "HeatMap", "HeatMap");
  26                 exit(1);
  27         }
  28 }
  29
  30 //**********************************************************************************************************************
  31
  32 string HeatMap::getPic(RAbundVector* rabund) {
  33         try {
  34
  35                 int numBinsToDisplay = rabund->getNumBins();
  36
  37                 if (numOTU != 0) { //user want to display a portion of the otus
  38                         if (numOTU < numBinsToDisplay) {  numBinsToDisplay = numOTU; }
  39                 }
  40
  41                 //sort lookup so shared bins are on top
  42                 if (sorted != "none") {  sortRabund(rabund);  }
  43
  44                 float maxRelAbund = 0.0;
  45
  46                 for(int i=0;i<rabund->size();i++){
  47                         float relAbund = rabund->get(i) / (float)rabund->getNumSeqs();
  48                         if(relAbund > maxRelAbund){     maxRelAbund = relAbund; }
  49                 }
  50
  51                 vector<string> scaleRelAbund(numBinsToDisplay, "");
  52
  53                 for(int i=0;i<numBinsToDisplay;i++){
  54                         float relAbund = rabund->get(i) / (float)rabund->getNumSeqs();
  55
  56                         if (m->control_pressed) { return "control"; }
  57
  58                         if (rabund->get(i) != 0) { //don't want log value of 0.
  59                                 if (scaler == "log10") {
  60                                         scaleRelAbund[i] = toHex(int(255 * log10(relAbund) / log10(maxRelAbund))) + "0000";
  61                                 }else if (scaler == "log2") {
  62                                         scaleRelAbund[i] = toHex(int(255 * log2(relAbund) / log2(maxRelAbund))) + "0000";
  63                                 }else if (scaler == "linear") {
  64                                         scaleRelAbund[i] = toHex(int(255 * relAbund / maxRelAbund)) + "0000";
  65                                 }else {  //if user enters invalid scaler option.
  66                                         scaleRelAbund[i] = toHex(int(255 * log10(relAbund / log10(maxRelAbund))))  + "0000";
  67                                 }
  68                         }
  69                         else { scaleRelAbund[i] = "FFFFFF";  }
  70                 }
  71
  72
  73                 string filenamesvg = outputDir + m->getRootName(m->getSimpleName(inputfile)) + rabund->getLabel() + ".heatmap.bin.svg";
  74                 m->openOutputFile(filenamesvg, outsvg);
  75
  76                 //svg image
  77                 outsvg << "<svg xmlns:svg=\"http://www.w3.org/2000/svg\" xmlns=\"http://www.w3.org/2000/svg\" width=\"100%\" height=\"100%\" viewBox=\"0 0 " + toString(300) + " " + toString((numBinsToDisplay*5 + 120))  + "\">\n";
  78                 outsvg << "<g>\n";
  79
  80                 //white backround
  81                 outsvg << "<rect fill=\"white\" stroke=\"white\" x=\"0\" y=\"0\" width=\"" + toString(300) + "\" height=\"" + toString((numBinsToDisplay*5 + 120))  + "\"/>";
  82                 outsvg << "<text fill=\"black\" class=\"seri\" text-anchor=\"middle\" font-size=\"" + toString(fontSize) + "\" x=\"" + toString((150) - 40) + "\" y=\"25\">Heatmap at distance " + rabund->getLabel() + "</text>\n";
  83
  84                 //output legend and color labels
  85                 string color;
  86                 int x = 0;
  87                 int y = 103 + (numBinsToDisplay*5);
  88                 printLegend(y, maxRelAbund);
  89
  90                 y = 70;
  91
  92                 for (int i = 0; i < scaleRelAbund.size(); i++) {
  93                         if (m->control_pressed) { outsvg.close(); return "control"; }
  94
  95                         outsvg << "<rect fill=\"#" + scaleRelAbund[i] + "\" stroke=\"#" + scaleRelAbund[i] + "\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"300\" height=\"5\"/>\n";
  96                         y += 5;
  97                 }
  98
  99                 outsvg << "</g>\n</svg>\n";
 100                 outsvg.close();
 101
 102                 return filenamesvg;
 103         }
 104         catch(exception& e) {
 105                 m->errorOut(e, "HeatMap", "getPic");
 106                 exit(1);
 107         }
 108 }
 109
 110 //**********************************************************************************************************************
 111
 112 string HeatMap::getPic(vector<SharedRAbundVector*> lookup) {
 113         try {
 114
 115                 int numBinsToDisplay = lookup[0]->size();
 116
 117                 if (numOTU != 0) { //user want to display a portion of the otus
 118                         if (numOTU < numBinsToDisplay) {  numBinsToDisplay = numOTU; }
 119                 }
 120
 121                 //sort lookup so shared bins are on top
 122                 if (sorted != "none") {  sortSharedVectors(lookup);  }
 123
 124                 vector<vector<string> > scaleRelAbund;
 125                 vector<float> maxRelAbund(lookup[0]->size(), 0.0);
 126                 float superMaxRelAbund = 0;
 127
 128                 for(int i = 0; i < lookup.size(); i++){
 129                         for(int j=0; j<lookup[i]->size(); j++){
 130
 131                                 float relAbund = lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs();
 132                                 if(relAbund > maxRelAbund[i]){  maxRelAbund[i] = relAbund;      }
 133                         }
 134                         if(maxRelAbund[i] > superMaxRelAbund){  superMaxRelAbund = maxRelAbund[i];      }
 135                 }
 136
 137                 scaleRelAbund.resize(lookup.size());
 138                 for(int i=0;i<lookup.size();i++){
 139                         scaleRelAbund[i].assign(numBinsToDisplay, "");
 140                         for(int j=0;j<numBinsToDisplay;j++){
 141                                 if (m->control_pressed) {  return "control"; }
 142                                 float relAbund = lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs();
 143
 144                                 if (lookup[i]->getAbundance(j) != 0) { //don't want log value of 0.
 145                                         if (scaler == "log10") {
 146                                                 scaleRelAbund[i][j] = toHex(int(255 * log10(relAbund) / log10(maxRelAbund[i]))) + "0000";
 147                                         }else if (scaler == "log2") {
 148                                                 scaleRelAbund[i][j] = toHex(int(255 * log2(relAbund) / log2(maxRelAbund[i]))) + "0000";
 149                                         }else if (scaler == "linear") {
 150                                                 scaleRelAbund[i][j] = toHex(int(255 * relAbund / maxRelAbund[i])) + "0000";
 151                                         }else {  //if user enters invalid scaler option.
 152                                                 scaleRelAbund[i][j] = toHex(int(255 * log10(relAbund / log10(maxRelAbund[i]))))  + "0000";
 153                                         }
 154                                 }else { scaleRelAbund[i][j] = "FFFFFF";  }
 155
 156                         }
 157                 }
 158
 159                 string filenamesvg = outputDir + m->getRootName(m->getSimpleName(inputfile)) + lookup[0]->getLabel() + ".heatmap.bin.svg";
 160                 m->openOutputFile(filenamesvg, outsvg);
 161
 162                 //svg image
 163                 outsvg << "<svg xmlns:svg=\"http://www.w3.org/2000/svg\" xmlns=\"http://www.w3.org/2000/svg\" width=\"100%\" height=\"100%\" viewBox=\"0 0 " + toString(lookup.size() * 300) + " " + toString((numBinsToDisplay*5 + 120))  + "\">\n";
 164                 outsvg << "<g>\n";
 165
 166                 //white backround
 167                 outsvg << "<rect fill=\"white\" stroke=\"white\" x=\"0\" y=\"0\" width=\"" + toString(lookup.size() * 300) + "\" height=\"" + toString((numBinsToDisplay*5 + 120))  + "\"/>";
 168                 outsvg << "<text fill=\"black\" class=\"seri\" font-size=\"" + toString(fontSize) + "\" text-anchor=\"middle\" x=\"" + toString((lookup.size() * 150) - 40) + "\" y=\"25\">Heatmap at distance " + lookup[0]->getLabel() + "</text>\n";
 169
 170                 //column labels
 171                 for (int h = 0; h < lookup.size(); h++) {
 172                         outsvg << "<text fill=\"black\" class=\"seri\" font-size=\"" + toString(fontSize) + "\" x=\"" + toString(((300 * (h+1)) - 150) - ((int)lookup[h]->getGroup().length() / 2)) + "\" y=\"50\">" + lookup[h]->getGroup() + "</text>\n";
 173                 }
 174
 175                 //output legend and color labels
 176                 string color;
 177                 int x = 0;
 178                 int y = 103 + (numBinsToDisplay*5);
 179                 printLegend(y, superMaxRelAbund);
 180
 181                 y = 70;
 182                 for (int i = 0; i < numBinsToDisplay; i++) {
 183                         for (int j = 0; j < scaleRelAbund.size(); j++) {
 184                                 if (m->control_pressed) { outsvg.close(); return "control"; }
 185
 186                                 outsvg << "<rect fill=\"#" + scaleRelAbund[j][i] + "\" stroke=\"#" + scaleRelAbund[j][i] + "\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"300\" height=\"5\"/>\n";
 187                                 x += 300;
 188                         }
 189                         x = 0;
 190                         y += 5;
 191                 }
 192
 193                 outsvg << "</g>\n</svg>\n";
 194                 outsvg.close();
 195
 196                 return filenamesvg;
 197
 198         }
 199         catch(exception& e) {
 200                 m->errorOut(e, "HeatMap", "getPic");
 201                 exit(1);
 202         }
 203 }
 204
 205 //**********************************************************************************************************************
 206 int HeatMap::sortSharedVectors(vector<SharedRAbundVector*>& lookup){
 207         try {
 208
 209                 vector<SharedRAbundVector*> looktemp;
 210                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 211                 map<int, int>::iterator it;
 212
 213                 /****************** find order of otus **********************/
 214                 if (sorted == "shared") {
 215                         place = orderShared(lookup);
 216                 }else if (sorted == "topotu") {
 217                         place = orderTopOtu(lookup);
 218                 }else if (sorted == "topgroup") {
 219                         place = orderTopGroup(lookup);
 220                 }else { m->mothurOut("Error: invalid sort option."); m->mothurOutEndLine();  return 1; }
 221
 222
 223                 /******************* create copy of lookup *********************/
 224                 //create and initialize looktemp as a copy of lookup
 225                 for (int i = 0; i < lookup.size(); i++) {
 226                         SharedRAbundVector* temp = new SharedRAbundVector(lookup[i]->getNumBins());
 227                         temp->setLabel(lookup[i]->getLabel());
 228                         temp->setGroup(lookup[i]->getGroup());
 229                         //copy lookup i's info
 230                         for (int j = 0; j < lookup[i]->size(); j++) {
 231                                 temp->set(j, lookup[i]->getAbundance(j), lookup[i]->getGroup());
 232                         }
 233                         looktemp.push_back(temp);
 234                 }
 235
 236                 /************************ fill lookup in order given by place *********************/
 237                 //for each bin
 238                 for (int i = 0; i < looktemp[0]->size(); i++) {                                                                                                         //place
 239                         //fill lookup                                                                                                                                                                   // 2 -> 1
 240                         for (int j = 0; j < looktemp.size(); j++) {                                                                                                             // 3 -> 2
 241                                 int newAbund = looktemp[j]->getAbundance(i);                                                                                            // 1 -> 3
 242                                 lookup[j]->set(place[i], newAbund, looktemp[j]->getGroup()); //binNumber, abundance, group
 243                         }
 244                 }
 245
 246                 //delete looktemp -- Sarah look at - this is causing segmentation faults
 247                 for (int j = 0; j < looktemp.size(); j++) {
 248 //                      delete looktemp[j];
 249                 }
 250
 251                 return 0;
 252
 253         }
 254         catch(exception& e) {
 255                 m->errorOut(e, "HeatMap", "sortSharedVectors");
 256                 exit(1);
 257         }
 258 }
 259 //**********************************************************************************************************************
 260 map<int, int> HeatMap::orderShared(vector<SharedRAbundVector*>& lookup){
 261         try {
 262
 263                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 264                 map<int, int>::iterator it;
 265
 266                 vector<int> sharedBins;
 267                 vector<int> uniqueBins;
 268
 269                 //for each bin
 270                 for (int i = 0; i < lookup[0]->size(); i++) {
 271                         int count = 0;
 272
 273                         //is this bin shared
 274                         for (int j = 0; j < lookup.size(); j++) {               if (lookup[j]->getAbundance(i) != 0) { count++; }       }
 275
 276                         if (count < 2)  {  uniqueBins.push_back(i); }
 277                         else                    {  sharedBins.push_back(i); }
 278                 }
 279
 280                 //fill place
 281                 for (int i = 0; i < sharedBins.size(); i++) { place[sharedBins[i]] = i; }
 282                 for (int i = 0; i < uniqueBins.size(); i++) { place[uniqueBins[i]] = (sharedBins.size() + i); }
 283
 284                 return place;
 285
 286         }
 287         catch(exception& e) {
 288                 m->errorOut(e, "HeatMap", "orderShared");
 289                 exit(1);
 290         }
 291 }
 292 //**********************************************************************************************************************
 293 map<int, int> HeatMap::orderTopOtu(vector<SharedRAbundVector*>& lookup){
 294         try {
 295
 296                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 297                 map<int, int>::iterator it;
 298
 299                 vector<binCount> totals;
 300
 301                 //for each bin
 302                 for (int i = 0; i < lookup[0]->size(); i++) {
 303                         int total = 0;
 304
 305                         for (int j = 0; j < lookup.size(); j++) {       total += lookup[j]->getAbundance(i);    }
 306
 307                         binCount temp(i, total);
 308
 309                         totals.push_back(temp);
 310                 }
 311
 312                 sort(totals.begin(), totals.end(), comparebinCounts);
 313
 314                 //fill place
 315                 for (int i = 0; i < totals.size(); i++) {   place[totals[i].bin] = i;  }
 316
 317                 return place;
 318
 319         }
 320         catch(exception& e) {
 321                 m->errorOut(e, "HeatMap", "orderTopOtu");
 322                 exit(1);
 323         }
 324 }
 325 //**********************************************************************************************************************
 326 map<int, int> HeatMap::orderTopGroup(vector<SharedRAbundVector*>& lookup){
 327         try {
 328
 329                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 330                 map<int, int>::iterator it;
 331
 332                 vector < vector<binCount> > totals; //totals[0] = bin totals for group 0, totals[1] = bin totals for group 1, ...
 333                 totals.resize(lookup.size());
 334
 335                 //for each bin
 336                 for (int i = 0; i < lookup[0]->size(); i++) {
 337                         for (int j = 0; j < lookup.size(); j++) {
 338                                 binCount temp(i, (lookup[j]->getAbundance(i)));
 339                                 totals[j].push_back(temp);
 340                         }
 341                 }
 342
 343                 for (int i = 0; i < totals.size(); i++) { sort(totals[i].begin(), totals[i].end(), comparebinCounts);  }
 344
 345                 //fill place
 346                 //grab the top otu for each group adding it if its not already added
 347                 int count = 0;
 348                 for (int i = 0; i < totals[0].size(); i++) {
 349
 350                         for (int j = 0; j < totals.size(); j++) {
 351                                 it = place.find(totals[j][i].bin);
 352
 353                                 if (it == place.end()) { //not added yet
 354                                         place[totals[j][i].bin] = count;
 355                                         count++;
 356                                 }
 357                         }
 358                 }
 359
 360                 return place;
 361
 362         }
 363         catch(exception& e) {
 364                 m->errorOut(e, "HeatMap", "orderTopGroup");
 365                 exit(1);
 366         }
 367 }
 368 //**********************************************************************************************************************
 369
 370 void HeatMap::printLegend(int y, float maxbin) {
 371         try {
 372
 373                 //output legend and color labels
 374                 //go through map and give each score a color value
 375                 string color;
 376                 int x = 10;
 377
 378                 //prints legend
 379                 for (int i = 1; i < 255; i++) {
 380                         color = toHex(int((float)(i)));
 381                         outsvg << "<rect fill=\"#" + color + "0000\" stroke=\"#" + color + "0000\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"1\" height=\"10\"/>\n";
 382                         x += 1;
 383                 }
 384
 385                 //prints legend labels
 386                 x = 10;
 387                 for (int i = 1; i<=5; i++) {
 388                         float label;
 389                         if(scaler== "log10")            {       label = maxbin * log10(51*i) / log10(255);      }
 390                         else if(scaler== "log2")        {       label = maxbin * log2(51*i) / log2(255);        }
 391                         else if(scaler== "linear")      {       label = maxbin * 51 * i / 255;                          }
 392                         else                                            {       label = maxbin * log10(51*i) / log10(255);      }
 393                         label = int(label * 1000 + 0.5);
 394                         label /= 1000.0;
 395                         string text = toString(label, 3);
 396
 397                         outsvg << "<text fill=\"black\" class=\"seri\" x=\"" + toString(x) + "\" y=\"" + toString(y-3) + "\">" + text + "</text>\n";
 398                         x += 60;
 399                 }
 400         }
 401
 402         catch(exception& e) {
 403                 m->errorOut(e, "HeatMap", "printLegend");
 404                 exit(1);
 405         }
 406 }
 407 //**********************************************************************************************************************
 408
 409 string HeatMap::getPic(vector<SharedRAbundFloatVector*> lookup) {
 410         try {
 411
 412                 int numBinsToDisplay = lookup[0]->size();
 413
 414                 if (numOTU != 0) { //user want to display a portion of the otus
 415                         if (numOTU < numBinsToDisplay) {  numBinsToDisplay = numOTU; }
 416                 }
 417
 418                 //sort lookup so shared bins are on top
 419                 if (sorted != "none") {  sortSharedVectors(lookup);  }
 420
 421                 vector<vector<string> > scaleRelAbund;
 422                 vector<float> maxRelAbund(lookup.size(), 0.0);
 423                 float superMaxRelAbund = 0;
 424
 425                 for(int i = 0; i < lookup.size(); i++){
 426                         for(int j=0; j<lookup[i]->size(); j++){
 427
 428                                 float relAbund = lookup[i]->getAbundance(j);
 429                                 if(relAbund > maxRelAbund[i]){  maxRelAbund[i] = relAbund;      }
 430                         }
 431                         if(maxRelAbund[i] > superMaxRelAbund){  superMaxRelAbund = maxRelAbund[i];      }
 432                 }
 433
 434                 scaleRelAbund.resize(lookup.size());
 435                 for(int i=0;i<lookup.size();i++){
 436                         scaleRelAbund[i].assign(numBinsToDisplay, "");
 437                         for(int j=0;j<numBinsToDisplay;j++){
 438                                 if (m->control_pressed) {  return "control"; }
 439                                 float relAbund = lookup[i]->getAbundance(j);
 440
 441                                 if (lookup[i]->getAbundance(j) != 0) { //don't want log value of 0.
 442                                         if (scaler == "log10") {
 443                                                 scaleRelAbund[i][j] = toHex(int(255 * log10(relAbund) / log10(maxRelAbund[i]))) + "0000";
 444                                         }else if (scaler == "log2") {
 445                                                 scaleRelAbund[i][j] = toHex(int(255 * log2(relAbund) / log2(maxRelAbund[i]))) + "0000";
 446                                         }else if (scaler == "linear") {
 447                                                 scaleRelAbund[i][j] = toHex(int(255 * relAbund / maxRelAbund[i])) + "0000";
 448                                         }else {  //if user enters invalid scaler option.
 449                                                 scaleRelAbund[i][j] = toHex(int(255 * log10(relAbund / log10(maxRelAbund[i]))))  + "0000";
 450                                         }
 451                                 }else { scaleRelAbund[i][j] = "FFFFFF";  }
 452
 453                         }
 454                 }
 455
 456                 string filenamesvg = outputDir + m->getRootName(m->getSimpleName(inputfile)) + lookup[0]->getLabel() + ".heatmap.bin.svg";
 457                 m->openOutputFile(filenamesvg, outsvg);
 458
 459                 //svg image
 460                 outsvg << "<svg xmlns:svg=\"http://www.w3.org/2000/svg\" xmlns=\"http://www.w3.org/2000/svg\" width=\"100%\" height=\"100%\" viewBox=\"0 0 " + toString(lookup.size() * 300) + " " + toString((numBinsToDisplay*5 + 120))  + "\">\n";
 461                 outsvg << "<g>\n";
 462
 463                 //white backround
 464                 outsvg << "<rect fill=\"white\" stroke=\"white\" x=\"0\" y=\"0\" width=\"" + toString(lookup.size() * 300) + "\" height=\"" + toString((numBinsToDisplay*5 + 120))  + "\"/>";
 465                 outsvg << "<text fill=\"black\" class=\"seri\" font-size=\"" + toString(fontSize) + "\" text-anchor=\"middle\" x=\"" + toString((lookup.size() * 150) - 40) + "\" y=\"25\">Heatmap at distance " + lookup[0]->getLabel() + "</text>\n";
 466
 467                 //column labels
 468                 for (int h = 0; h < lookup.size(); h++) {
 469                         outsvg << "<text fill=\"black\" class=\"seri\" font-size=\"" + toString(fontSize) + "\" x=\"" + toString(((300 * (h+1)) - 150) - ((int)lookup[h]->getGroup().length() / 2)) + "\" y=\"50\">" + lookup[h]->getGroup() + "</text>\n";
 470                 }
 471
 472                 //output legend and color labels
 473                 string color;
 474                 int x = 0;
 475                 int y = 103 + (numBinsToDisplay*5);
 476                 printLegend(y, superMaxRelAbund);
 477
 478                 y = 70;
 479                 for (int i = 0; i < numBinsToDisplay; i++) {
 480                         for (int j = 0; j < scaleRelAbund.size(); j++) {
 481                                 if (m->control_pressed) { outsvg.close(); return "control"; }
 482
 483                                 outsvg << "<rect fill=\"#" + scaleRelAbund[j][i] + "\" stroke=\"#" + scaleRelAbund[j][i] + "\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"300\" height=\"5\"/>\n";
 484                                 x += 300;
 485                         }
 486                         x = 0;
 487                         y += 5;
 488                 }
 489
 490                 outsvg << "</g>\n</svg>\n";
 491                 outsvg.close();
 492
 493                 return filenamesvg;
 494
 495         }
 496         catch(exception& e) {
 497                 m->errorOut(e, "HeatMap", "getPic");
 498                 exit(1);
 499         }
 500 }
 501 //**********************************************************************************************************************
 502 int HeatMap::sortSharedVectors(vector<SharedRAbundFloatVector*>& lookup){
 503         try {
 504
 505                 vector<SharedRAbundFloatVector*> looktemp;
 506                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 507                 map<int, int>::iterator it;
 508
 509                 /****************** find order of otus **********************/
 510                 if (sorted == "shared") {
 511                         place = orderShared(lookup);
 512                 }else if (sorted == "topotu") {
 513                         place = orderTopOtu(lookup);
 514                 }else if (sorted == "topgroup") {
 515                         place = orderTopGroup(lookup);
 516                 }else { m->mothurOut("Error: invalid sort option."); m->mothurOutEndLine();  return 1; }
 517
 518
 519                 /******************* create copy of lookup *********************/
 520                 //create and initialize looktemp as a copy of lookup
 521                 for (int i = 0; i < lookup.size(); i++) {
 522                         SharedRAbundFloatVector* temp = new SharedRAbundFloatVector(lookup[i]->getNumBins());
 523                         temp->setLabel(lookup[i]->getLabel());
 524                         temp->setGroup(lookup[i]->getGroup());
 525                         //copy lookup i's info
 526                         for (int j = 0; j < lookup[i]->size(); j++) {
 527                                 temp->set(j, lookup[i]->getAbundance(j), lookup[i]->getGroup());
 528                         }
 529                         looktemp.push_back(temp);
 530                 }
 531
 532                 /************************ fill lookup in order given by place *********************/
 533                 //for each bin
 534                 for (int i = 0; i < looktemp[0]->size(); i++) {                                                                                                         //place
 535                         //fill lookup                                                                                                                                                                   // 2 -> 1
 536                         for (int j = 0; j < looktemp.size(); j++) {                                                                                                             // 3 -> 2
 537                                 float newAbund = looktemp[j]->getAbundance(i);                                                                                          // 1 -> 3
 538                                 lookup[j]->set(place[i], newAbund, looktemp[j]->getGroup()); //binNumber, abundance, group
 539                         }
 540                 }
 541
 542                 //delete looktemp -- Sarah look at - this is causing segmentation faults
 543                 for (int j = 0; j < looktemp.size(); j++) {
 544 //                      delete looktemp[j];
 545                 }
 546
 547                 return 0;
 548
 549         }
 550         catch(exception& e) {
 551                 m->errorOut(e, "HeatMap", "sortSharedVectors");
 552                 exit(1);
 553         }
 554 }
 555 //**********************************************************************************************************************
 556 int HeatMap::sortRabund(RAbundVector*& r){
 557         try {
 558                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 559                 map<int, int>::iterator it;
 560
 561                 /****************** find order of otus **********************/
 562                 vector<binCount> totals;
 563
 564                 //for each bin
 565                 for (int i = 0; i < r->getNumBins(); i++) {
 566                         binCount temp(i, r->get(i));
 567
 568                         totals.push_back(temp);
 569                 }
 570
 571                 sort(totals.begin(), totals.end(), comparebinCounts);
 572
 573                 //fill place
 574                 for (int i = 0; i < totals.size(); i++) {   place[totals[i].bin] = i;  }
 575
 576                 /******************* create copy of lookup *********************/
 577                 //create and initialize rtemp as a copy of r
 578
 579                 RAbundVector* rtemp = new RAbundVector(r->getNumBins());
 580                 for (int i = 0; i < r->size(); i++) {  rtemp->set(i, r->get(i)); }
 581                 rtemp->setLabel(r->getLabel());
 582
 583                 /************************ fill lookup in order given by place *********************/
 584                 //for each bin
 585                 for (int i = 0; i < rtemp->size(); i++) {                                                                               //place
 586                         //fill lookup                                                                                                                           // 2 -> 1
 587                                                                                                                                                                                 // 3 -> 2
 588                         int newAbund = rtemp->get(i);                                                                                           // 1 -> 3
 589                         r->set(place[i], newAbund); //binNumber, abundance
 590                 }
 591
 592                 return 0;
 593
 594         }
 595         catch(exception& e) {
 596                 m->errorOut(e, "HeatMap", "sortRabund");
 597                 exit(1);
 598         }
 599 }
 600 //**********************************************************************************************************************
 601 map<int, int> HeatMap::orderShared(vector<SharedRAbundFloatVector*>& lookup){
 602         try {
 603
 604                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 605                 map<int, int>::iterator it;
 606
 607                 vector<int> sharedBins;
 608                 vector<int> uniqueBins;
 609
 610                 //for each bin
 611                 for (int i = 0; i < lookup[0]->size(); i++) {
 612                         int count = 0;
 613
 614                         //is this bin shared
 615                         for (int j = 0; j < lookup.size(); j++) {               if (lookup[j]->getAbundance(i) != 0) { count++; }       }
 616
 617                         if (count < 2)  {  uniqueBins.push_back(i); }
 618                         else                    {  sharedBins.push_back(i); }
 619                 }
 620
 621                 //fill place
 622                 for (int i = 0; i < sharedBins.size(); i++) { place[sharedBins[i]] = i; }
 623                 for (int i = 0; i < uniqueBins.size(); i++) { place[uniqueBins[i]] = (sharedBins.size() + i); }
 624
 625                 return place;
 626
 627         }
 628         catch(exception& e) {
 629                 m->errorOut(e, "HeatMap", "orderShared");
 630                 exit(1);
 631         }
 632 }
 633 //**********************************************************************************************************************
 634 map<int, int> HeatMap::orderTopOtu(vector<SharedRAbundFloatVector*>& lookup){
 635         try {
 636
 637                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 638                 map<int, int>::iterator it;
 639
 640                 vector<binCountFloat> totals;
 641
 642                 //for each bin
 643                 for (int i = 0; i < lookup[0]->size(); i++) {
 644                         int total = 0;
 645
 646                         for (int j = 0; j < lookup.size(); j++) {       total += lookup[j]->getAbundance(i);    }
 647
 648                         binCountFloat temp(i, total);
 649
 650                         totals.push_back(temp);
 651                 }
 652
 653                 sort(totals.begin(), totals.end(), comparebinFloatCounts);
 654
 655                 //fill place
 656                 for (int i = 0; i < totals.size(); i++) {   place[totals[i].bin] = i;  }
 657
 658                 return place;
 659
 660         }
 661         catch(exception& e) {
 662                 m->errorOut(e, "HeatMap", "orderTopOtu");
 663                 exit(1);
 664         }
 665 }
 666 //**********************************************************************************************************************
 667 map<int, int> HeatMap::orderTopGroup(vector<SharedRAbundFloatVector*>& lookup){
 668         try {
 669
 670                 map<int, int> place; //spot in lookup where you insert shared by, ie, 3 -> 2 if they are shared by 3 inset into location 2.
 671                 map<int, int>::iterator it;
 672
 673                 vector < vector<binCountFloat> > totals; //totals[0] = bin totals for group 0, totals[1] = bin totals for group 1, ...
 674                 totals.resize(lookup.size());
 675
 676                 //for each bin
 677                 for (int i = 0; i < lookup[0]->size(); i++) {
 678                         for (int j = 0; j < lookup.size(); j++) {
 679                                 binCountFloat temp(i, (lookup[j]->getAbundance(i)));
 680                                 totals[j].push_back(temp);
 681                         }
 682                 }
 683
 684                 for (int i = 0; i < totals.size(); i++) { sort(totals[i].begin(), totals[i].end(), comparebinFloatCounts);  }
 685
 686                 //fill place
 687                 //grab the top otu for each group adding it if its not already added
 688                 int count = 0;
 689                 for (int i = 0; i < totals[0].size(); i++) {
 690
 691                         for (int j = 0; j < totals.size(); j++) {
 692                                 it = place.find(totals[j][i].bin);
 693
 694                                 if (it == place.end()) { //not added yet
 695                                         place[totals[j][i].bin] = count;
 696                                         count++;
 697                                 }
 698                         }
 699                 }
 700
 701                 return place;
 702
 703         }
 704         catch(exception& e) {
 705                 m->errorOut(e, "HeatMap", "orderTopGroup");
 706                 exit(1);
 707         }
 708 }
 709 //**********************************************************************************************************************
 710
 711
 712
 713
 714