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(){
  14         try {
  15                 globaldata = GlobalData::getInstance();
  16                 format = globaldata->getFormat();
  17                 sorted = globaldata->getSorted();
  18                 util = new SharedUtil();
  19
  20         }
  21         catch(exception& e) {
  22                 cout << "Standard Error: " << e.what() << " has occurred in the HeatMap class Function HeatMap. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
  23                 exit(1);
  24         }
  25         catch(...) {
  26                 cout << "An unknown error has occurred in the HeatMap class function HeatMap. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
  27                 exit(1);
  28         }
  29 }
  30 //**********************************************************************************************************************
  31 void HeatMap::getPic(OrderVector* order) {
  32         try {
  33                 colorScale.clear();
  34
  35                 rabund = order->getRAbundVector();
  36
  37                 //get users scaling method
  38                 scaler = globaldata->getScale();
  39
  40                 float maxbin = 0.0;
  41                 for (int i = 0; i < rabund.size(); i++) {
  42                         if (rabund.get(i) != 0) { //don't want log value of 0.
  43                                         if (scaler == "log10") {
  44                                                 colorScale[(log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))] = "";
  45                                                 if (maxbin < (log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))) { maxbin = (log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000)); }
  46                                         }else if (scaler == "log2") {
  47                                                 colorScale[(log2((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))] = "";
  48                                                 if (maxbin < (log2((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))) { maxbin = (log2((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000)); }
  49                                         }else if (scaler == "linear") {
  50                                                 colorScale[rabund.get(i)] = "";
  51                                                 if (maxbin < rabund.get(i)) { maxbin = rabund.get(i); }
  52                                         }else {  //if user enters invalid scaler option.
  53                                                 cout << scaler << " is not a valid scaler option. I will use log10." << endl;
  54                                                 colorScale[(log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))] = "";
  55                                                 if (maxbin < (log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))) { maxbin = (log10((rabund.get(i)) / (float)rabund.getNumSeqs()) * 1000); }
  56                                         }
  57                         }else { colorScale[0] = "00";  }
  58                 }
  59
  60                 float scalers = 255 / (float) maxbin;
  61
  62                 //go through map and give each score a color value
  63                 for (it = colorScale.begin(); it != colorScale.end(); it++) {
  64                         it->second = toHex(int(float(it->first) * scalers));
  65                         if(it->second.length() == 1) {  it->second = "0" + it->second;  }
  66                 }
  67
  68                 string filenamesvg = getRootName(globaldata->inputFileName) + order->getLabel() + ".heatmap.svg";
  69                 openOutputFile(filenamesvg, outsvg);
  70
  71                 //svg image
  72                 outsvg << "<svg width=\"100%\" height=\"100%\" viewBox=\"0 0 300 " + toString((rabund.getNumBins()*5 + 120))  + "\">\n";
  73                 outsvg << "<g>\n";
  74
  75                 //white backround
  76                 outsvg << "<rect fill=\"white\" stroke=\"white\" x=\"0\" y=\"0\" width=\"300\" height=\"" + toString((rabund.getNumBins()*5 + 120))  + "\"/>";
  77                 outsvg << "<text fill=\"black\" class=\"seri\" x=\"100\" y=\"25\">Heatmap at distance " + order->getLabel() + "</text>\n";
  78
  79                 //output legend and color labels
  80                 string color;
  81                 int x = 0;
  82                 int y = 103 + (rabund.getNumBins()*5);
  83                 printLegend(y, maxbin);
  84
  85                 y = 70;
  86                 for (int i = 0; i <= rabund.getNumBins(); i++) {
  87                         if (rabund.get(i) != 0) { //don't want log value of 0.
  88                                 if (scaler == "log10") {
  89                                         color = colorScale[(log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))];
  90                                 }else if (scaler == "log2") {
  91                                         color = colorScale[(log2((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))];
  92                                 }else if (scaler == "linear") {
  93                                         color = colorScale[rabund.get(i)];
  94                                 }else {  color = colorScale[(log10((rabund.get(i) / (float)rabund.getNumSeqs()) * 1000))];      }
  95                         }else { color = "OO";  }
  96
  97                         outsvg << "<rect fill=\"#" + color + "0000\" stroke=\"#" + color + "0000\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"300\" height=\"5\"/>\n";
  98                         y += 5;
  99                 }
 100
 101                 outsvg << "</g>\n</svg>\n";
 102                 outsvg.close();
 103
 104         }
 105         catch(exception& e) {
 106                 cout << "Standard Error: " << e.what() << " has occurred in the HeatMap class Function getPic. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 107                 exit(1);
 108         }
 109         catch(...) {
 110                 cout << "An unknown error has occurred in the HeatMap class function getPic. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 111                 exit(1);
 112         }
 113 }
 114 //**********************************************************************************************************************
 115 void HeatMap::getPic(SharedOrderVector* sharedorder) {
 116         try {
 117                 colorScale.clear();
 118
 119                 //fills vector of sharedsabunds - lookup
 120                 util->getSharedVectors(globaldata->Groups, lookup, sharedorder);  //fills group vectors from order vector.
 121
 122                 //sort lookup so shared bins are on top
 123                 if (sorted == "T") {  sortSharedVectors();  }
 124
 125                 //get users scaling method
 126                 scaler = globaldata->getScale();
 127
 128                 float maxbin = 0.0;
 129                 for (int i = 0; i < lookup.size(); i++) {
 130                         for (int j = 0; j < lookup[i]->size(); j++) {
 131                                 if (lookup[i]->getAbundance(j) != 0) { //don't want log value of 0.
 132                                         if (scaler == "log10") {
 133                                                 colorScale[(log10((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000))] = "";
 134                                                 if (maxbin < (log10((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000))) { maxbin = (log10((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000)); }
 135                                         }else if (scaler == "log2") {
 136                                                 colorScale[(log2((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000))] = "";
 137                                                 if (maxbin < (log2((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000))) { maxbin = (log2((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000)); }
 138                                         }else if (scaler == "linear") {
 139                                                 colorScale[lookup[i]->getAbundance(j)] = "";
 140                                                 if (maxbin < lookup[i]->getAbundance(j)) { maxbin = lookup[i]->getAbundance(j); }
 141                                         }else {  //if user enters invalid scaler option.
 142                                                 cout << scaler << " is not a valid scaler option. I will use log10." << endl;
 143                                                 colorScale[(log10((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000))] = "";
 144                                                 if (maxbin < (log10((lookup[i]->getAbundance(j) / (float)lookup[i]->getNumSeqs()) * 1000))) { maxbin = (log10((lookup[i]->getAbundance(j)) / (float)lookup[i]->getNumSeqs()) * 1000); }
 145                                         }
 146                                 }else { colorScale[0] = "00";  }
 147                         }
 148                 }
 149
 150                 //get scaler
 151                 float scalers = 255 / (float) maxbin;
 152
 153
 154                 //go through map and give each score a color value
 155                 for (it = colorScale.begin(); it != colorScale.end(); it++) {
 156                         it->second = toHex(int(float(it->first) * scalers));
 157                         if(it->second.length() == 1) {  it->second = "0" + it->second;  }
 158                 }
 159
 160                 string filenamesvg = getRootName(globaldata->inputFileName) + sharedorder->getLabel() + ".heatmap.svg";
 161                 openOutputFile(filenamesvg, outsvg);
 162
 163                 //svg image
 164                 outsvg << "<svg width=\"100%\" height=\"100%\" viewBox=\"0 0 " + toString(lookup.size() * 300) + " " + toString((lookup[0]->getNumBins()*5 + 120))  + "\">\n";
 165                 outsvg << "<g>\n";
 166
 167                 //white backround
 168                 outsvg << "<rect fill=\"white\" stroke=\"white\" x=\"0\" y=\"0\" width=\"" + toString(lookup.size() * 300) + "\" height=\"" + toString((lookup[0]->getNumBins()*5 + 120))  + "\"/>";
 169                 outsvg << "<text fill=\"black\" class=\"seri\" x=\"" + toString((lookup.size() * 150) - 40) + "\" y=\"25\">Heatmap at distance " + sharedorder->getLabel() + "</text>\n";
 170
 171                 //column labels
 172                 for (int h = 0; h < lookup.size(); h++) {
 173                         outsvg << "<text fill=\"black\" class=\"seri\" x=\"" + toString(((300 * (h+1)) - 150) - ((int)lookup[h]->getGroup().length() / 2)) + "\" y=\"50\">" + lookup[h]->getGroup() + "</text>\n";
 174                 }
 175
 176                 //output legend and color labels
 177                 string color;
 178                 int x = 0;
 179                 int y = 103 + (lookup[0]->getNumBins()*5);
 180                 printLegend(y, maxbin);
 181
 182                 y = 70;
 183                 for (int i = 0; i < lookup[0]->size(); i++) {
 184                         for (int j = 0; j < lookup.size(); j++) {
 185
 186                                 if (lookup[j]->getAbundance(i) != 0) { //don't want log value of 0.
 187                                         if (scaler == "log10") {
 188                                                 color = colorScale[(log10((lookup[j]->getAbundance(i) / (float)lookup[j]->getNumSeqs()) * 1000))];
 189                                         }else if (scaler == "log2") {
 190                                                 color = colorScale[(log2((lookup[j]->getAbundance(i) / (float)lookup[j]->getNumSeqs()) * 1000))];
 191                                         }else if (scaler == "linear") {
 192                                                 color = colorScale[lookup[j]->getAbundance(i)];
 193                                         }else {  color = colorScale[(log10((lookup[j]->getAbundance(i) / (float)lookup[j]->getNumSeqs()) * 1000))];     }
 194                                 }else { color = "OO";  }
 195
 196
 197                                 outsvg << "<rect fill=\"#" + color + "0000\" stroke=\"#" + color + "0000\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"300\" height=\"5\"/>\n";
 198                                 x += 300;
 199                         }
 200                         x = 0;
 201                         y += 5;
 202                 }
 203
 204
 205                 outsvg << "</g>\n</svg>\n";
 206                 outsvg.close();
 207
 208         }
 209         catch(exception& e) {
 210                 cout << "Standard Error: " << e.what() << " has occurred in the HeatMap class Function getPic. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 211                 exit(1);
 212         }
 213         catch(...) {
 214                 cout << "An unknown error has occurred in the HeatMap class function getPic. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 215                 exit(1);
 216         }
 217 }
 218
 219 //**********************************************************************************************************************
 220 void HeatMap::sortSharedVectors(){
 221         try {
 222                 //copy lookup and then clear it to refill with sorted.
 223                 //loop though lookup and determine if they are shared
 224                 //if they are then insert in the front
 225                 //if not push to back
 226
 227                 vector<SharedRAbundVector*> looktemp;
 228                 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.
 229                 map<int, int>::iterator it;
 230                 int count;
 231
 232                 //create and initialize looktemp as a copy of lookup
 233                 for (int i = 0; i < lookup.size(); i++) {
 234                         SharedRAbundVector* temp = new SharedRAbundVector(lookup[i]->getNumBins());
 235                         temp->setLabel(lookup[i]->getLabel());
 236                         temp->setGroup(lookup[i]->getGroup());
 237                         //copy lookup i's info
 238                         for (int j = 0; j < lookup[i]->size(); j++) {
 239                                 temp->set(j, lookup[i]->getAbundance(j), lookup[i]->getGroup());
 240                         }
 241                         looktemp.push_back(temp);
 242                 }
 243
 244                 //clear out lookup to create sorted lookup
 245                 lookup.clear();
 246
 247                 //create and initialize lookup to empty vectors
 248                 for (int i = 0; i < looktemp.size(); i++) {
 249                         SharedRAbundVector* temp = new SharedRAbundVector();
 250                         temp->setLabel(looktemp[i]->getLabel());
 251                         temp->setGroup(looktemp[i]->getGroup());
 252                         lookup.push_back(temp);
 253
 254                         //initialize place map
 255                         place[i] = 0;
 256                 }
 257
 258
 259                 //for each bin
 260                 for (int i = 0; i < looktemp[0]->size(); i++) {
 261                         count = 0;
 262                         bool updatePlace = false;
 263                         //for each group
 264                         for (int j = 0; j < looktemp.size(); j++) {
 265                                 if (looktemp[j]->getAbundance(i) != 0) { count++; }
 266                         }
 267
 268                         //fill lookup
 269                         for (int j = 0; j < looktemp.size(); j++) {
 270                                 //if they are not shared then push to back, if they are not insert in front
 271                                 if (count < 2)  { lookup[j]->push_back(looktemp[j]->getAbundance(i), i, looktemp[j]->getGroup()); }
 272                                 //they are shared by some
 273                                 else {  lookup[j]->insert(looktemp[j]->getAbundance(i), place[count], looktemp[j]->getGroup());   updatePlace = true; }
 274                         }
 275
 276                         if (updatePlace == true) {
 277                                 //move place holders below where you entered up to "make space" for you entry
 278                                 for (it = place.begin(); it!= place.end(); it++) {
 279                                         if (it->first < count) { it->second++; }
 280                                 }
 281                         }
 282                 }
 283
 284                 //delete looktemp
 285                 for (int j = 0; j < looktemp.size(); j++) {
 286                         delete looktemp[j];
 287                 }
 288
 289         }
 290         catch(exception& e) {
 291                 cout << "Standard Error: " << e.what() << " has occurred in the HeatMap class Function sortSharedVectors. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 292                 exit(1);
 293         }
 294         catch(...) {
 295                 cout << "An unknown error has occurred in the HeatMap class function sortSharedVectors. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 296                 exit(1);
 297         }
 298
 299 }
 300
 301 //**********************************************************************************************************************
 302 void HeatMap::printLegend(int y, float maxbin) {
 303         try {
 304
 305                 //output legend and color labels
 306                 //go through map and give each score a color value
 307                 string color;
 308                 int x = 0;
 309                 if (maxbin != 0) {
 310                         //convert maxbin to relative abundance again
 311                         if (scaler == "log10") {
 312                                 maxbin = pow(10, maxbin) / 1000;
 313                         }else if (scaler == "log2") {
 314                                 maxbin = pow(2, maxbin) / 1000;
 315                         }else {  maxbin = pow(10, maxbin) / 1000;       }
 316                 }else { maxbin = 0.00; }
 317
 318                 //5 is the number of boxes in the legend
 319                 float maxbinScaler = maxbin / 10;
 320                 float colorScaler = 255 / 10;
 321
 322                 //prints legend
 323                 for (int i = 0; i < 10; i++) {
 324                         color = toHex(int((float)(i+1) * colorScaler));
 325                         outsvg << "<rect fill=\"#" + color + "0000\" stroke=\"#" + color + "0000\" x=\"" + toString(x) + "\" y=\"" + toString(y) + "\" width=\"30\" height=\"10\"/>\n";
 326                         x += 30;
 327                 }
 328
 329                 //prints legend labels
 330                 x -= 30;
 331                 for (int i = 10; i > 0; i = i-2) {
 332                         string label = toString((i * maxbinScaler));
 333                         //set precision of relative abundance to 3
 334                         int pos = label.find_first_of('.');
 335                         label = label.substr(0,pos+4);
 336
 337                         outsvg << "<text fill=\"black\" class=\"seri\" x=\"" + toString(x) + "\" y=\"" + toString(y-3) + "\">" + label + "</text>\n";
 338                         x -= 60;
 339                 }
 340         }
 341         catch(exception& e) {
 342                 cout << "Standard Error: " << e.what() << " has occurred in the HeatMap class Function printLegend. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 343                 exit(1);
 344         }
 345         catch(...) {
 346                 cout << "An unknown error has occurred in the HeatMap class function printLegend. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n";
 347                 exit(1);
 348         }
 349
 350 }
 351
 352
 353
 354