src/hg/instinct/hgHeatmap2/drawingCode.c 1.66
1.66 2009/06/04 03:50:36 jsanborn
added copyright notices, removed cluster library
Index: src/hg/instinct/hgHeatmap2/drawingCode.c
===================================================================
RCS file: /projects/compbio/cvsroot/kent/src/hg/instinct/hgHeatmap2/drawingCode.c,v
retrieving revision 1.65
retrieving revision 1.66
diff -b -B -U 1000000 -r1.65 -r1.66
--- src/hg/instinct/hgHeatmap2/drawingCode.c 6 May 2009 00:06:08 -0000 1.65
+++ src/hg/instinct/hgHeatmap2/drawingCode.c 4 Jun 2009 03:50:36 -0000 1.66
@@ -1,3299 +1,3302 @@
+/********************************************************************************/
+/* Copyright 2007-2009 -- The Regents of the University of California */
+/********************************************************************************/
+
/* mainPage - drs the main hgHeatmap page, including some controls
* on the top and the graphic. */
#define EXPR_DATA_SHADES 16
#define DEFAULT_MAX_DEVIATION 0.7
#define COLOR_SCALE 1
#define MAX_PCA_COMPONENTS 2
//#define DBL_MAX 1000000.0 // Changed for 32-bit computers compatibility
#define PROBE_MARK_HEIGHT 3
#define PROBE_MARK_BUFFER 2
#include <limits.h>
#include "common.h"
#include "bed.h"
#include "cart.h"
#include "customTrack.h"
#include "errCatch.h"
#include "genoLay.h"
#include "trackLayout.h"
#include "hash.h"
#include "hdb.h"
#include "hgColors.h"
#include "hPrint.h"
#include "htmshell.h"
#include "jsHelper.h"
#include "psGfx.h"
#include "trashDir.h"
#include "vGfx.h"
#include "hvGfx.h"
#include "web.h"
#include "cytoBand.h"
#include "hCytoBand.h"
#include "hgHeatmapLib.h"
#include "heatmapUtility.h"
#include "featuresLib.h"
#include "hgStatsLib.h"
#include "hgHeatmap2.h"
#include "filterFeatures.h"
#include "chromGraph.h"
#include "hgChromGraph.h"
-#include "hgStats.h"
static char const rcsid[] = "$Id$";
static char *heatMapDbProfile = "localDb"; // database profile to use
struct gifList {
struct gifList *next;
int expId;
char *filename;
};
struct chromZoom {
struct chromZoom *next;
char *fullName;
int size;
int baseStart;
int baseEnd;
};
struct trackLayout tl; /* Dimensions of things, fonts, etc. */
char *blockStatTest[3]; /* compatability with library code */
#define CLIP(p,limit) if (p < 0) p = 0; if (p >= (limit)) p = (limit)-1;
void verticalTextCentered(struct vGfx *vg, int x, int y,
int width, int height, int colorIx, MgFont *font, char *string)
/* Draw a vertical line of text in middle of the box.
* The string will read from bottom to top
* This is not perfect by any means, but seems to draw the labels fairly nice
* for the pathways.
*/
{
/* don't let this run wild */
CLIP(width, vg->width);
CLIP(height, vg->height);
if ((width <= 0) || (height <= 0))
return;
struct vGfx *vgHoriz;
int i, j;
/* reversed meanings of width and height here since this is going
* to rotate 90 degrees
*/
int offset = 0;
int cutoff = height / 7; // TODO, approximately correct.
int lineSpacing = 15;
vgHoriz = vgOpenGif(height, width, "/dev/null");
struct slName *sl, *slList = slNameListFromString(string, '_');
struct dyString *dy = newDyString(0), *dyList = NULL;
for (sl = slList; sl; sl = sl->next)
{ /* loop through words in text, separated by underscores */
if ( (dy->stringSize > 0) && (offset < width) &&
(dy->stringSize + strlen(sl->name) > cutoff))
{ /* enough string for a single line is ready to be drawn and can be drawn */
slAddHead(&dyList, dy);
dy = newDyString(0);
offset += lineSpacing;
}
dyStringAppend(dy, sl->name);
dyStringAppend(dy, " ");
}
if ( (dy->stringSize > 0) && (offset < width) )
/* Clean up remaining string */
slAddHead(&dyList, dy);
slReverse(&dyList);
int count = slCount(dyList);
offset = 0;
for (dy = dyList; dy ; dy = dy->next)
{ /* -3 is needed to approximately recenter text */
vgTextCentered(vgHoriz, 0, -3*(count - 1), height, width+offset,
colorIx, font, dy->string);
offset += lineSpacing;
}
dyStringFreeList(&dyList);
/* now, blit from the horizontal to the vertical, rotate -90 (CCW) */
for (i = 0; i < height; ++i) /* xSrc -> yDest */
{
int yDest = height - i;
for (j = 0; j < width; ++j) /* ySrc -> xDest */
vgDot(vg, j+x, yDest+y, vgGetDot(vgHoriz, i, j));
}
vgClose(&vgHoriz);
}
void verticalTextRight(struct vGfx *vg, int x, int y,
int width, int height, int colorIx, MgFont *font, char *string)
/* Draw a vertical line of text in middle of the box.
* The string will read from bottom to top
* This is not perfect by any means, but seems to draw the labels fairly nice
* for the pathways.
*/
{
/* don't let this run wild */
CLIP(width, vg->width);
CLIP(height, vg->height);
if ((width <= 0) || (height <= 0))
return;
struct vGfx *vgHoriz;
int i, j;
/* reversed meanings of width and height here since this is going
* to rotate 90 degrees
*/
int offset = 0;
int cutoff = height / 7; // TODO, approximately correct.
int lineSpacing = 15;
vgHoriz = vgOpenGif(height, width, "/dev/null");
string = replaceChars(string, "_", " ");
struct slName *sl, *slList = slNameListFromString(string, ' ');
struct dyString *dy = newDyString(0), *dyList = NULL;
for (sl = slList; sl; sl = sl->next)
{ /* loop through words in text, separated by underscores */
if ( (dy->stringSize > 0) && (offset < width) &&
(dy->stringSize + strlen(sl->name) > cutoff))
{ /* enough string for a single line is ready to be drawn and can be drawn */
slAddHead(&dyList, dy);
dy = newDyString(0);
offset += lineSpacing;
}
dyStringAppend(dy, sl->name);
dyStringAppend(dy, " ");
}
if ( (dy->stringSize > 0) && (offset < width) )
/* Clean up remaining string */
slAddHead(&dyList, dy);
slReverse(&dyList);
int count = slCount(dyList);
offset = 0;
for (dy = dyList; dy ; dy = dy->next)
{ /* -3 is needed to approximately recenter text */
vgTextRight(vgHoriz, 0, -3*(count - 1), height + 2, width+offset,
colorIx, font, dy->string);
offset += lineSpacing;
}
dyStringFreeList(&dyList);
/* now, blit from the horizontal to the vertical, rotate -90 (CCW) */
for (i = 0; i < height; ++i) /* xSrc -> yDest */
{
int yDest = height - i;
for (j = 0; j < width; ++j) /* ySrc -> xDest */
vgDot(vg, j+x, yDest+y, vgGetDot(vgHoriz, i, j));
}
vgClose(&vgHoriz);
}
boolean singleChromImage(struct heatmapLay *hl)
{
if (!hl)
return FALSE;
if (slCount(hl->elements) == 1)
return TRUE;
return FALSE;
}
/* return an array for reordering the experiments in a chromosome */
double maxDeviation(char* heatmap)
{
struct hashEl *e = hashLookup(ghHash, heatmap);
if (!e)
errAbort("Could not found heatmap %s", heatmap);
struct genoHeatmap *gh = e->val;
if (gh->expScale >0)
return gh->expScale;
return DEFAULT_MAX_DEVIATION;
}
void vgMakeColorGradient(struct vGfx *vg,
struct rgbColor *start, struct rgbColor *end,
int steps, Color *colorIxs)
/* Make a color gradient that goes smoothly from start
* to end colors in given number of steps. Put indicesgl->chromLis
* in color table in colorIxs */
{
double scale = 0, invScale;
double invStep;
int i; int r,g,b;
steps -= 1; /* Easier to do the calculation in an inclusive way. */
invStep = 1.0/steps;
for (i=0; i<=steps; ++i)
{
invScale = 1.0 - scale;
r = invScale * start->r + scale * end->r;
g = invScale * start->g + scale * end->g;
b = invScale * start->b + scale * end->b;
colorIxs[i] = vgFindColorIx(vg, r, g, b);
scale += invStep;
}
}
char *getId(struct sqlConnection *conn, char *table, char *key, char *sample, char *value)
/* get patient ID from sample (or experiment) Id */
{
char query[512];
safef(query, sizeof(query), "select %s from %s where %s = '%s' ", key, table, value, sample);
return sqlQuickString(conn, query);
}
char *getPathwayDb()
/* Get pathway database */
{
return "pathway";
}
struct sqlConnection *getPathwayDbConn()
/* Get connection to pathway database */
{
struct sqlConnection *conn = hAllocConnProfile(heatMapDbProfile, getPathwayDb());
return conn;
}
struct sqlConnection *getFeatureDbConn(struct genoHeatmap *gh)
/* Get feature database connection */
{
char *featureDb = gh->patDb;
char *raName = gh->raFile;
if ((featureDb) && (raName))
return hAllocConnProfile(heatMapDbProfile, featureDb);
else
return NULL;
}
void heatmapLaySetColor(struct heatmapLay *hlList, struct rgbColor *low, struct rgbColor *zero,
struct rgbColor *high)
{
struct heatmapLay *hl;
for (hl = hlList; hl; hl = hl->next)
{
hl->lowColor = low;
hl->zeroColor = zero;
hl->highColor = high;
}
}
void heatmapLaySetDefaultColor(struct heatmapLay *hl, char *platform)
{
struct rgbColor *low = NULL;
struct rgbColor *zero = NULL;
struct rgbColor *high = NULL;
// default, expression data color scheme:
if (sameString(platform, "CNV") || sameString(platform, "SNP"))
{ // CNV, SNP data
AllocVar(zero); // white
zero->r = 255;
zero->g = 255;
zero->b = 255;
AllocVar(high); // red
high->r = 255;
high->g = 0;
high->b = 0;
AllocVar(low); // blue
low->r = 0;
low->g = 0;
low->b = 255;
}
else
{ // expression data, others...
AllocVar(zero); // black
zero->r = 0;
zero->g = 0;
zero->b = 0;
AllocVar(high); // red
high->r = 255;
high->g = 0;
high->b = 0;
AllocVar(low); // green
low->r = 0;
low->g = 255;
low->b = 0;
}
heatmapLaySetColor(hl, low, zero, high);
}
void heatmapLayResetMinMax(struct heatmapLay *hlList)
{
if (!hlList)
return;
struct heatmapLay *hl;
for (hl = hlList; hl; hl = hl->next)
{
hl->minVal = INT_MAX;
hl->minValPos = -1;
hl->maxVal = INT_MIN;
hl->maxValPos = -1;
}
}
void drawLayoutLines(struct vGfx *vg, struct heatmapLay *hlList, int width)
{
if (!hlList)
return;
struct heatmapLay *hl;
for (hl = hlList; hl; hl = hl->next)
{
int leftX = hl->pixelStart - 1;
if (leftX > 0)
vgBox(vg, leftX, 0, 1, hlList->height, MG_GRAY);
}
}
struct heatmapLay *chromLayout(struct sqlConnection *conn, struct genoHeatmap *gh,
char *chromStart, int startBase,
char *chromEnd, int endBase,
int width, int sampleHeight)
{
/* Get list of chromosomes and lay them out. */
struct sqlConnection *hg18conn = hAllocConnProfile(heatMapDbProfile, database);
struct genoLayChrom *chrom, *allChromList = genoLayDbChroms(hg18conn, FALSE);
if (!allChromList)
return NULL;
unsigned long totalBases = 0;
int chromCount = 0;
int keepChroms = 0;
if (!chromStart)
chromStart = cloneString(allChromList->fullName);
if (!chromEnd)
for (chrom = allChromList; chrom; chrom = chrom->next)
if (!chrom->next)
chromEnd = cloneString(chrom->fullName);
struct chromZoom *cz, *czList = NULL;
for (chrom = allChromList; chrom != NULL; chrom = chrom->next)
{
if (sameString(chrom->fullName, chromStart))
keepChroms = 1;
if (keepChroms)
{
AllocVar(cz);
cz->fullName = cloneString(chrom->fullName);
cz->baseStart = -1;
cz->baseEnd = -1;
if (sameString(chrom->fullName, chromStart) && (startBase != -1))
cz->baseStart = startBase;
if (sameString(chrom->fullName, chromEnd) && (endBase != -1))
{
if (endBase > chrom->size)
endBase = chrom->size;
cz->baseEnd = endBase + 1;
}
if (cz->baseStart < 0)
cz->baseStart = 0;
if (cz->baseEnd < 0)
cz->baseEnd = chrom->size;
cz->size = cz->baseEnd - cz->baseStart;
slAddHead(&czList, cz);
totalBases += cz->size;
chromCount += 1;
}
if (sameString(chrom->fullName, chromEnd))
keepChroms = 0;
}
slReverse(&czList);
int availablePixels = width - (chromCount - 1);
double basesPerPixel = (double) totalBases / availablePixels;
struct heatmapLay *hm = AllocA(struct heatmapLay);
struct hmElement *hEl, *heList = NULL;
hm->name = cloneString("chromInfo");
hm->height = 0;
hm->sampleHeight = sampleHeight;
hm->width = width;
hm->basesPerPixel = basesPerPixel;
trackLayoutInit(&tl, cart);
hm->fontHeight = tl.fontHeight;
hm->font = tl.font;
hm->minVal = INT_MAX;
hm->minValPos = -1;
hm->maxVal = INT_MIN;
hm->maxValPos = -1;
int xOff = 0;
for (cz = czList; cz != NULL; cz = cz->next)
{
unsigned long bStart = cz->baseStart;
unsigned long bEnd = cz->baseEnd;
int pixStart = round(0 / basesPerPixel) + xOff;
int pixEnd = round((bEnd-bStart) / basesPerPixel) + xOff;
if (pixEnd >= width)
pixEnd = width;
AllocVar(hEl);
hEl->name = cloneString(cz->fullName);
hEl->pixelStart = pixStart;
hEl->pixelEnd = pixEnd;
hEl->baseStart = bStart;
hEl->baseEnd = bEnd;
hEl->shortLabel = NULL;
hEl->longLabel = NULL;
hEl->probeName = NULL;
slAddHead(&heList, hEl);
xOff = pixEnd + 1;
}
slReverse(&heList);
hm->elements = heList;
if (gh)
heatmapLaySetDefaultColor(hm, gh->platform);
return hm;
}
struct heatmapLay *genesetLayout(struct genoHeatmap *gh,
char *geneSetNames,
int width, int sampleHeight)
{
struct geneSet *geneSets = getAllPathways(cart, getPathwayDb(), geneSetNames);
geneSetsChromLayout(geneSets);
boolean isSNP = FALSE;
struct hash *geneHash = NULL;
if (gh)
{
getChromHeatmapHash(&geneHash, gh->database, gh->probeTable,
gh->name, NULL, geneSets);
if (!geneHash)
return NULL;
if (sameString(gh->platform, "SNP"))
isSNP = TRUE;
}
else
{ /* for geneset labels, gh == NULL, but still want layout) */
geneHash = hashNew(0);
}
geneSetsAddPathwayLayout(geneSets, geneHash);
trackLayoutInit(&tl, cart);
struct geneSet *gs;
struct heatmapLay *hl, *hlList = NULL;
for (gs = geneSets; gs; gs = gs->next)
{
AllocVar(hl);
hl->name = cloneString(gs->name);
hl->height = 0;
hl->sampleHeight = sampleHeight;
hl->width = gs->width;
hl->fontHeight = tl.fontHeight;
hl->font = tl.font;
hl->pixelStart = gs->x;
hl->pixelEnd = gs->x + gs->width;
hl->minVal = INT_MAX;
hl->minValPos = -1;
hl->maxVal = INT_MIN;
hl->maxValPos = -1;
struct slName *sl;
double count = 0.0;
int gsStart = hl->pixelStart;
struct hmElement *hEl, *heList = NULL;
for (sl = gs->genes; sl; sl = sl->next)
{
int numProbes = 0;
struct hashEl *el = hashLookup(geneHash, sl->name);
while (el)
{
el = hashLookupNext(el);
numProbes++;
}
int geneStart = (int) ( gs->pixelsPerGene * count + gsStart );
if (numProbes == 0)
{
int pixStart = geneStart;
int pixEnd = pixStart + ceil(gs->pixelsPerGene);
AllocVar(hEl);
hEl->name = cloneString(sl->name);
hEl->pixelStart = pixStart;
hEl->pixelEnd = pixEnd;
hEl->baseStart = 0;
hEl->baseEnd = 0;
hEl->shortLabel = NULL;
hEl->longLabel = NULL;
hEl->probeName = cloneString("N/A");
hEl->bed = NULL;
slAddHead(&heList, hEl);
}
else
{
el = hashLookup(geneHash, sl->name);
double probeCount = 0.0;
double pixelsPerProbe = gs->pixelsPerGene / (double) numProbes;
while (el)
{
struct bed *nb = el->val;
el = hashLookupNext(el);
int pixStart, pixEnd;
if (isSNP)
{
pixStart = geneStart;
pixEnd = pixStart + ceil(gs->pixelsPerGene);
}
else
{
pixStart = (int) ( pixelsPerProbe * probeCount + geneStart );
pixEnd = pixStart + ceil(pixelsPerProbe);
}
AllocVar(hEl);
hEl->name = cloneString(sl->name);
hEl->pixelStart = pixStart;
hEl->pixelEnd = pixEnd;
hEl->baseStart = nb->chromStart;
hEl->baseEnd = nb->chromEnd;
hEl->shortLabel = NULL;
hEl->longLabel = NULL;
hEl->probeName = cloneString(nb->name);
hEl->bed = nb;
slAddHead(&heList, hEl);
probeCount += 1.0;
}
}
count += 1.0;
}
slReverse(&heList);
hl->elements = heList;
slAddHead(&hlList, hl);
}
slReverse(&hlList);
if (gh)
heatmapLaySetDefaultColor(hlList, gh->platform);
return hlList;
}
struct heatmapLay *featureLayout(struct genoHeatmap *gh, int width, int sampleHeight)
/* Layout Feature Sorter based on previously layed out heatmaps */
{
struct sqlConnection *conn = getFeatureDbConn(gh);
if (!conn)
return NULL;
char *raName = gh->raFile;
struct column *col, *colList = getColumns(conn, raName, gh->patDb);
hFreeConn(&conn);
int numVis = 0;
for (col = colList; col != NULL; col = col->next)
if (col->on)
numVis++;
int availablePixels = width;
double basesPerPixel = (double) numVis / availablePixels;
struct heatmapLay *hm = AllocA(struct heatmapLay);
struct hmElement *hEl, *heList = NULL;
hm->name = cloneString("featureInfo");
hm->height = 0;
hm->sampleHeight = sampleHeight;
hm->width = width;
hm->basesPerPixel = basesPerPixel;
trackLayoutInit(&tl, cart);
hm->font = tl.font;
hm->fontHeight = tl.fontHeight;
double pixelsPerBase = 0.0;
if (basesPerPixel > 0)
pixelsPerBase = 1.0/basesPerPixel;
int colIx = 0;
for (col = colList; col != NULL; col = col->next)
{
if (!col->on)
continue;
AllocVar(hEl);
hEl->name = cloneString(col->name);
hEl->shortLabel = cloneString(col->shortLabel);
hEl->longLabel = cloneString(col->longLabel);
hEl->pixelStart = round(pixelsPerBase * colIx);
hEl->pixelEnd = round(pixelsPerBase * (colIx + 1));
hEl->baseStart = -1;
hEl->baseEnd = -1;
hEl->probeName = NULL;
colIx++;
slAddHead(&heList, hEl);
}
slReverse(&heList);
hm->elements = heList;
return hm;
}
void geneSetsAddPathwayLayout(struct geneSet *geneSets, struct hash *geneHash)
{
/* lay out of all pathways, add setting the gs->pixelsPerGene gs->activeNum for each gene */
if (geneSets == NULL || geneHash == NULL)
return;
struct geneSet *gs;
for (gs = geneSets; gs; gs = gs->next)
{
gs->numGenesActive = slCount(gs->genes);//numActive;
gs->pixelsPerGene = (double) gs->width / (double) gs->numGenesActive;
}
}
void geneSetsChromLayout(struct geneSet *geneSets)
/* lay out of all pathways similar to chromsomes,
without setting the gs->pixelsPerGene gs->activeNum for each gene */
{
if (geneSets == NULL)
return;
int totalGenes = 0;
struct geneSet *gs;
int numSets = 0;
for (gs = geneSets; gs; gs = gs->next)
{
totalGenes += gs->numGenes;
numSets++;
}
if (totalGenes == 0)
return;
int picWidth = cartUsualInt(cart, hghImageWidth, hgHeatmapDefaultPixWidth);
double xOff = 0;
int minWidth = 20;
int totalW = (double) (picWidth - numSets);
for (gs = geneSets; gs; gs = gs->next)
{
int w = gs->numGenes * (picWidth - numSets) / totalGenes;
if (w < minWidth)
{
totalGenes -= gs->numGenes;
totalW -= minWidth;
}
}
/* 5 is hard-coded , it works well*/
double pixelsPerGene = (totalW - 5 ) / (double) totalGenes;
for (gs = geneSets; gs; gs = gs->next)
{
gs->x = (int) xOff;
gs->y = 0;
gs->width = (int) (gs->numGenes * pixelsPerGene);
if (gs->width < minWidth)
gs->width = minWidth;
xOff += (double) gs->width + 1.0;
}
}
void drawGeneSetSingleStat(struct vGfx *vg, struct hash *statHash, struct heatmapLay *hl,
double gMin, double gMax, float colorCutoff,
char *chromHeatmap)
{
if (!statHash | !hl)
return;
Color valCol;
double height = hghBed5Height;
double baseline = (float) height /2.0;
double gScale = baseline;
if (gMax != gMin)
gScale /= (gMax - gMin);
vgBox(vg, hl->pixelStart, baseline-1, hl->width, 1, MG_GRAY);
int leftX = hl->pixelStart - 1;
if (leftX < 0)
leftX = 0;
int rightX = hl->pixelStart + hl->width;
vgBox(vg, leftX, 0, 1, height, MG_GRAY);
vgBox(vg, rightX, 0, 1, height, MG_GRAY);
float val, dirVal;
int x;
int y1,y2;
int w, h;
struct hmElement *hEl;
for (hEl = hl->elements; hEl ; hEl = hEl->next)
{
struct hgStats *nb = NULL;
struct hashEl *el;
el = hashLookup(statHash, hEl->name);
while (el)
{ // TODO: do this in a smarter way.
nb = el->val;
if (sameString(nb->name, hEl->probeName))
break;
nb = NULL;
el = hashLookupNext(el);
}
if (!nb)
continue;
nb = el->val;
x = hEl->pixelStart; // get the X position
w = hEl->pixelEnd - hEl->pixelStart;
val = fabs(nb->prob);
dirVal = nb->stats;
if (val>gMax)
val = gMax;
if (val <gMin)
val = gMin;
if (val < colorCutoff)
valCol = MG_GRAY;
else if (dirVal > 0)
valCol = MG_RED;
else
valCol = MG_GREEN;
double tmpV;
double tmpH;
if (dirVal>0)
{
y1 = baseline - val *gScale; //smaller value
y2 = baseline; //greater value
h = y2-y1;
/* sometimes if the bed segment is too small, resulting in width=0,
* reset to 1 pixel length */
if ((w ==0) && (h != 0))
w=1;
if (w && h)
vgBox(vg,x,y1, w, h,valCol);
tmpV = val;
tmpH = y1;
}
else
{
y1 = baseline + val *gScale; //greater value
y2 = baseline; //smaller value
h = y1-y2;
/* sometimes if the bed segment is too small, resulting in width=0,
* reset to 1 pixel length */
if ((w ==0) && (h != 0))
w=1;
if (w && h)
vgBox(vg,x,y2,w,h,valCol);
tmpV = -1.0 * val;
tmpH = y2 + h;
}
/* For scaling */
if (tmpV > hl->maxVal)
{
hl->maxVal = tmpV;
hl->maxValPos = tmpH;
}
if (tmpV < hl->minVal)
{
hl->minVal = tmpV;
hl->minValPos = tmpH;
}
el = hashLookupNext(el);
}
}
/* Draw hgStats results along gene sets */
void drawGeneSetStats(struct vGfx *vg, struct hash *statHash,
struct heatmapLay *hlList, double gMin, double gMax, float colorCutoff,
char *tableName)
/* Drawing code */
{
if (!statHash)
errAbort("No hgStats results");
struct heatmapLay *hl;
for (hl = hlList; hl ; hl = hl->next)
drawGeneSetSingleStat(vg, statHash, hl, gMin, gMax, colorCutoff, tableName);
}
void drawChromStats(struct vGfx *vg, struct hgStats *bg, struct heatmapLay *hl,
double gMin, double gMax, double colorCutoff)
{
if (!bg)
errAbort("No hgStats data");
double height = hghBed5Height;
double baseline = (float) height / 2.0;
vgBox(vg, 0, baseline-1, hl->width, 1, MG_GRAY);
/* Chromosome layout hashs */
struct hmElement *hEl = NULL;
struct hash *chromHashX = newHash(0);
struct hash *pixelHashX = newHash(0);
for(hEl = hl->elements; hEl; hEl = hEl->next)
{
int chromX = hEl->baseStart;
int pixelX = hEl->pixelStart;
char *name = hEl->name;
hashAddInt(chromHashX, name, chromX);
hashAddInt(pixelHashX, name, pixelX);
int leftX = hEl->pixelStart - 1;
if (leftX < 0)
leftX = 0;
int rightX = hEl->pixelEnd;
if (rightX > hl->width - 1)
rightX = hl->width - 1;
vgBox(vg, leftX, 0, 1, height, MG_GRAY);
vgBox(vg, rightX, 0, 1, height, MG_GRAY);
}
double pixelsPerBase = 1.0/hl->basesPerPixel;
Color valCol;
double gScale = height / 2.0;
if (gMax != gMin)
gScale /= (gMax - gMin);
int start,end;
char *chromName;
float val, dirVal;
int x1,x2;
int y1;
int w, h;
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hmPixel *hm = NULL;
struct hgStats *ibg = NULL;
for(ibg = bg; ibg ; ibg = ibg->next)
{
chromName = ibg->chrom;
start = ibg->chromStart;
end = ibg->chromEnd;
if (!hashLookup(chromHashX, chromName))
continue;
int chromX = hashIntVal(chromHashX, chromName);
int pixelX = hashIntVal(pixelHashX, chromName);
x1 = pixelsPerBase*(start - chromX) + pixelX;
x2 = pixelsPerBase*(end - chromX) + pixelX;
w = x2-x1;
val = fabs(ibg->prob); // probability score control the absolute height of the drawing code
if (val > gMax)
val = gMax;
if (val < gMin)
val = gMin;
dirVal = ibg->stats; // statisitcs control drawing above or below baseline and color
if (dirVal > 0)
y1 = baseline - val *gScale; //smaller value
else
y1 = baseline; //smaller value
h = val*gScale;
double tmpV = val;
double tmpH = y1;
if (dirVal < 0)
{
tmpV = -1.0 * tmpV;
tmpH = tmpH + h;
}
if (tmpV > hl->maxVal)
{
hl->maxVal = tmpV;
hl->maxValPos = tmpH;
}
if (tmpV < hl->minVal)
{
hl->minVal = tmpV;
hl->minValPos = tmpH;
}
if (dirVal > 0)
safef(pixelStr,sizeof (pixelStr),"%d,%d+",x1,x2);
else
safef(pixelStr,sizeof (pixelStr),"%d,%d-",x1,x2);
struct hashEl *el = hashLookup(pixelHash, pixelStr);
if (!el)
{
hm = AllocA(struct hmPixel);
hm->x = x1;
hm->y = y1;
hm->w = w;
hm->h = h;
hm->val = dirVal;
hashAdd(pixelHash, pixelStr, hm);
}
else
{
hm = el->val;
//TO DOS: the current implementation is taking the most significant value,
// a more sensible way to do it probably should be taking the average.
// Implement in the future.
if (hm->h <h)
{
hm->y = y1;
hm->w = w;
hm->h = h;
hm->val = dirVal;
}
}
}
struct hashEl *elList = hashElListHash(pixelHash);
struct hashEl *el;
for (el = elList; el; el=el->next)
{
hm = el->val;
w = hm->w;
h = hm->h;
x1 = hm->x;
y1 = hm->y;
val = hm->val;
/* sometimes if the bed segment is too small, resulting in width=0,
* reset to 1 pixel length
*/
if ((w ==0) && (h != 0))
w = 1;
if (w < 1 && h < 1)
continue;
if ( h < colorCutoff*gScale)
valCol= MG_GRAY;
else if (val > 0)
valCol = MG_RED;
else
valCol = MG_GREEN;
vgBox(vg, x1, y1, w, h, valCol);
}
hashElFreeList(&elList);
}
void drawSubgroups(struct vGfx *vg, char *chromHeatmap, int width,
int height, int sampleHeight)
/* Draw features to right of chromosomes in layout at given
* y offset and height. */
{
if (!chromHeatmap)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
char *raName = gh->raFile;
int subsetNum = cartUsualInt(cart, hgh2SubgroupNum, hgh2SubgroupDefaultNum);
/* This is the key function to call for subgrouping */
struct slName **ptSubsets = getSubsets(gh, subsetNum, raName );
if (!ptSubsets)
return;
int subset = 0;
Color valCol;
char *sample;
double pixelsPerSubgroup = (double) width / (double) subsetNum;
for (subset = 0; subset < subsetNum; subset++)
{
if (!ptSubsets[subset])
continue;
struct slName *subSampleList = ptSubsets[subset];
int x = round(pixelsPerSubgroup * subset);
int w = ceil(pixelsPerSubgroup);
vgSetClip(vg, x, 0, w, height);
vgBox(vg, x, 0, w, height, MG_WHITE);
struct slName *sl = NULL;
for (sl = gh->sampleList; sl ; sl = sl->next)
{
sample = sl->name;
int orderId = hashIntValDefault(gh->sampleOrder, sample, -1);
if (orderId == -1)
continue;
if (slNameInList( subSampleList, sample))
valCol = advFilterColor (subset);
else
continue;
int h = sampleHeight;
int y = orderId * h;
vgBox(vg, x, y, w, h, valCol);
}
}
for (subset=0; subset < subsetNum; subset++)
if (ptSubsets[subset])
slNameFreeList(ptSubsets[subset]);
free(ptSubsets);
vgUnclip(vg);
}
int ifDrawFeatureLabel(struct heatmapLay *hl)
{
if (!hl)
return 0;
struct hmElement *hEl = hl->elements;
if (!hEl)
return 0;
int width = hEl->pixelEnd - hEl->pixelStart;
if (width < (hl->fontHeight + 5))
return 0;
// Feature width is wide enough to accomodate text.
return 1;
}
void drawFeatures(struct vGfx *vg, char *chromHeatmap, struct heatmapLay *hl, int yOff,
boolean reverse)
/* Draw features to right of chromosomes in layout at given
* y offset and height. */
{
if (!hl)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
/* set up connection to table with both patient and sample information */
char *labTable = gh->patTable;
char *value = gh->sampleField;
char *key = gh->patField;
char *db = gh->patDb;
if ((labTable == NULL) || (key == NULL) || (value == NULL) || (db==NULL))
return;
struct sqlConnection *conn = hAllocConnProfile(heatMapDbProfile, db);
if (!conn)
return;
double colorScale = 0.0;
double val;
double absVal, minVal, maxVal;
char *minCutVal, *maxCutVal;
double offset;
int reverseColor;
char *colorScheme;
Color valCol;
Color upShades[EXPR_DATA_SHADES];
Color downShades[EXPR_DATA_SHADES];
static struct rgbColor black = {0, 0, 0};
static struct rgbColor red = {255, 0, 0};
static struct rgbColor green = {0, 255, 0};
static struct rgbColor yellow = {255, 255, 0};
struct slName *id = NULL;
char *raName = gh->raFile;
struct column *col, *colList = getColumns(conn, raName, gh->patDb);
int sampleHeight = hl->sampleHeight;
int height = heatmapHeight(gh) * sampleHeight;
struct hmElement *hEl;
for (hEl = hl->elements; hEl; hEl = hEl->next)
{
int chromX = hEl->pixelStart, chromY = 0;
int width = hEl->pixelEnd - hEl->pixelStart;
int h = sampleHeight;
for (col = colList; col != NULL; col = col->next)
if (sameString(col->name, hEl->name))
break;
if (!col)
continue;
colorScheme = col->cellColorSchemeVal(col);
if sameWord(colorScheme,"redGreen")
{
vgMakeColorGradient(vg, &black, &red, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, &black, &green, EXPR_DATA_SHADES, downShades);
}
if sameWord(colorScheme,"default")
{
vgMakeColorGradient(vg, &black, &yellow, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, &black, &green, EXPR_DATA_SHADES, downShades);
}
if (ifDrawFeatureLabel(hl))
verticalTextRight(vg, chromX, height, width, hghFeatureLabel,
MG_BLACK, hl->font, col->shortLabel);
vgSetClip(vg, chromX, chromY, width, height);
vgBox(vg, chromX, chromY, width, height, MG_GRAY);
offset = atof(col->cellOffsetVal(col));
reverseColor = atoi(col->cellColorReverseVal(col));
char *minValStr = col->cellMinVal(col, conn);
if (!minValStr)
minVal = 0.0;
else
minVal = atof(minValStr) + offset;
char *maxValStr = col->cellMaxVal(col, conn);
if (!maxValStr)
maxVal = 0.0;
else
maxVal = atof(maxValStr) + offset;
if (reverseColor == -1)
{
maxVal = reverseColor*maxVal;
minVal = reverseColor*maxVal;
}
minCutVal = col->cellMinCutVal(col, conn);
maxCutVal = col->cellMaxCutVal(col, conn);
/* double color: consistent with the scale of microarray data display */
if ((minVal<0) && (maxVal>0))
{
double md = max(fabs(minVal), maxVal);
colorScale = COLOR_SCALE / md;
}
else /* single color */
colorScale = COLOR_SCALE / (maxVal - minVal);
struct slName *sl = NULL;
for (sl = gh->sampleList; sl ; sl = sl->next)
{
int orderId = hashIntValDefault(gh->sampleOrder, sl->name, -1);
if (orderId == -1)
continue;
id = slNameNew(getId(conn, labTable, key, sl->name, value));
char *cellVal = col->cellVal(col, id, conn);
if (!cellVal)
continue;
valCol = MG_GRAY;
val = atof(cellVal);
if (minCutVal)
if (val < atof(minCutVal))
continue;
if (maxCutVal)
if (val > atof(maxCutVal))
continue;
absVal = fabs(val);
absVal = reverseColor *( absVal +offset);
/* double color: consistant with the scale of microarray data display */
int colorIndex ;
if ((minVal<0) && (maxVal>0))
colorIndex = (int)((absVal) * (EXPR_DATA_SHADES-1.0) * colorScale);
else
colorIndex = (int)((absVal-minVal) * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if (reverse)
val = val * -1.0;
if (val >= 0.0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
int y = chromY + orderId * h;
vgBox(vg, chromX, y, width, h, valCol);
freez(&id);
}
vgUnclip(vg);
}
hFreeConn(&conn);
}
int hmPixelCmpVal(const void *va, const void *vb)
/* Compare to sort columns based on priority. */
{
const struct hmPixel *a = *((struct hmPixel **)va);
const struct hmPixel *b = *((struct hmPixel **)vb);
float dif = a->val/a->count - b->val/b->count;
if (dif < 0)
return -1;
else if (dif > 0)
return 1;
else
return 0;
}
void drawFeatureSubgroupBar(struct vGfx *vg, struct slName *sampleList, int sampleHeight,
boolean isBottomSubgroup)
{
if (!sampleList)
return;
int height = hgHeatmapDefaultSummaryHeight;
int width = hghSubgroupDefaultPixWidth / 2;
vgSetClip(vg, 0, 0, hghSubgroupDefaultPixWidth, height);
if (isBottomSubgroup) // only draw labels on bottom/green subgroup
vgBox(vg, width, 0, width, height, MG_GREEN);
else
vgBox(vg, 0, 0, width, height, MG_RED);
vgUnclip(vg);
}
void drawFeatureSummary(struct vGfx *vg, char *chromHeatmap, struct heatmapLay *hl, int yOff,
boolean reverse, struct slName *sampleList, boolean isBottomSubgroup)
/* Draw feature summary to right of chromosomes in layout at given
* y offset and height... actually a feature heatmap where every column is sorted ascending
* breaking all relationships between horizontal rows
*/
{
if (!hl)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
/* set up connection to table with both patient and sample information */
char *labTable = gh->patTable;
char *value = gh->sampleField;
char *key = gh->patField;
char *db = gh->patDb;
struct slName *sl;
struct slInt *si, *siList = NULL;
for (sl = sampleList; sl; sl = sl->next)
{
int i = slNameFindIx(gh->sampleList, sl->name);
if (i < 0)
continue;
si = slIntNew(i);
slAddHead(&siList, si);
}
slReverse(&siList);
if ((labTable == NULL) || (key == NULL) || (value == NULL) || (db==NULL))
return;
struct sqlConnection *conn = hAllocConnProfile(heatMapDbProfile, db);
if (!conn)
return;
double colorScale = 0.0;
double val;
double absVal, minVal, maxVal;
char *minCutVal, *maxCutVal;
double offset;
int reverseColor;
char *colorScheme;
Color valCol;
Color upShades[EXPR_DATA_SHADES];
Color downShades[EXPR_DATA_SHADES];
static struct rgbColor black = {0, 0, 0};
static struct rgbColor red = {255, 0, 0};
static struct rgbColor green = {0, 255, 0};
static struct rgbColor yellow = {255, 255, 0};
struct slName *id = NULL;
char *raName = gh->raFile;
struct column *col, *colList = getColumns(conn, raName, gh->patDb);
int height = hl->height;
int buffer = 1;
double h = 1.0;
if (slCount(sampleList))
h = (double) hl->height / (double) slCount(sampleList);
struct hmElement *hEl;
for (hEl = hl->elements; hEl; hEl = hEl->next)
{
int chromX = hEl->pixelStart;
int width = hEl->pixelEnd - hEl->pixelStart;
for (col = colList; col != NULL; col = col->next)
if (sameString(col->name, hEl->name))
break;
if (!col)
continue;
colorScheme = col->cellColorSchemeVal(col);
if sameWord(colorScheme,"redGreen")
{
vgMakeColorGradient(vg, &black, &red, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, &black, &green, EXPR_DATA_SHADES, downShades);
}
if sameWord(colorScheme,"default")
{
vgMakeColorGradient(vg, &black, &yellow, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, &black, &green, EXPR_DATA_SHADES, downShades);
}
if (ifDrawFeatureLabel(hl) && isBottomSubgroup)
verticalTextRight(vg, chromX, height, width, hghFeatureLabel,
MG_BLACK, hl->font, col->shortLabel);
/* To differentiate between normal feature heatmaps and summary heatmaps
* keep a one-pixel white buffer between each feature column, done by
* reducing width of clipped area */
vgSetClip(vg, chromX-buffer, 0, width-buffer, height);
vgBox(vg, chromX, 0, width, height, MG_GRAY);
offset = atof(col->cellOffsetVal(col));
reverseColor = atoi(col->cellColorReverseVal(col));
char *minValStr = col->cellMinVal(col, conn);
if (!minValStr)
minVal = 0.0;
else
minVal = atof(minValStr) + offset;
char *maxValStr = col->cellMaxVal(col, conn);
if (!maxValStr)
maxVal = 0.0;
else
maxVal = atof(maxValStr) + offset;
if (reverseColor == -1)
{
maxVal = reverseColor*maxVal;
minVal = reverseColor*maxVal;
}
minCutVal = col->cellMinCutVal(col, conn);
maxCutVal = col->cellMaxCutVal(col, conn);
/* double color: consistent with the scale of microarray data display */
if ((minVal<0) && (maxVal>0))
{
double md = max(fabs(minVal), maxVal);
colorScale = COLOR_SCALE / md;
}
else /* single color */
colorScale = COLOR_SCALE / (maxVal - minVal);
struct slName *sl;
struct hmPixel *hm, *hmList = NULL;
for (sl = sampleList; sl ; sl = sl->next)
{
int orderId = hashIntValDefault(gh->sampleOrder, sl->name, -1);
if (orderId == -1)
continue;
int y = round((double) orderId * h);
id = slNameNew(getId(conn, labTable, key, sl->name, value));
char *cellVal = col->cellVal(col, id, conn);
if (!cellVal)
continue;
val = atof(cellVal);
if (minCutVal)
if (val < atof(minCutVal))
continue;
if (maxCutVal)
if (val > atof(maxCutVal))
continue;
absVal = fabs(val);
absVal = reverseColor *( absVal +offset);
/* double color: consistant with the scale of microarray data display */
int colorIndex ;
if ((minVal<0) && (maxVal>0))
colorIndex = (int)((absVal) * (EXPR_DATA_SHADES-1.0) * colorScale);
else
colorIndex = (int)((absVal-minVal) * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if (reverse)
val = val * -1.0;
if (val >= 0.0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
AllocVar(hm);
hm->val = val;
hm->count = 1;
hm->x = chromX;
hm->y = y;
hm->w = width;
hm->h = ceil(h);
hm->color = valCol;
slAddHead(&hmList, hm);
freez(&id);
}
/* Sort column by value */
slSort(&hmList, hmPixelCmpVal);
int yOff = 0;
/* Draw sorted data */
for (hm = hmList; hm; hm = hm->next)
{
vgBox(vg, hm->x, yOff, hm->w, hm->h, hm->color);
yOff += ceil(h);
}
vgUnclip(vg);
slFreeList(&hmList);
}
hFreeConn(&conn);
}
void drawGeneSetHeatmapsByPixel(struct vGfx *vg, char* database,
char *chromHeatmap, struct heatmapLay *hlList)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
if (!hlList)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
int *chromOrder = getBedOrder(gh);
struct bed *nb = NULL;
float gain = gh->gainSet;
double md = maxDeviation(chromHeatmap);
double colorScale = COLOR_SCALE / md;
double val;
double absVal;
int valId;
Color valCol;
Color upShades[EXPR_DATA_SHADES];
Color downShades[EXPR_DATA_SHADES];
vgMakeColorGradient(vg, hlList->zeroColor, hlList->highColor, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, hlList->zeroColor, hlList->lowColor, EXPR_DATA_SHADES, downShades);
int sampleHeight = hlList->sampleHeight;
int height = heatmapHeight(gh) * sampleHeight;
int maxWidth = 0;
struct heatmapLay *hl;
for (hl = hlList; hl ; hl = hl->next)
maxWidth = hl->pixelEnd;
Color zeroColor = vgFindColorIx(vg, hlList->zeroColor->r,
hlList->zeroColor->g, hlList->zeroColor->b);
boolean isSNP = FALSE;
if (sameString(gh->platform, "SNP"))
isSNP = TRUE;
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hmPixel *hm, *hmList = NULL;
for (hl = hlList; hl ; hl = hl->next)
{
int width = hl->pixelEnd - hl->pixelStart;
vgSetClip(vg, hl->pixelStart, 0, width, height);
vgBox(vg, hl->pixelStart, 0, maxWidth, height, zeroColor);
struct hmElement *hEl;
for (hEl = hl->elements; hEl ; hEl = hEl->next)
{
nb = hEl->bed;
if (!nb)
continue;
int x = hEl->pixelStart;
int w = hEl->pixelEnd - hEl->pixelStart;
int h = sampleHeight;
int i;
for(i = 0; i < nb->expCount; ++i)
{
val = nb->expScores[i];
valId = nb->expIds[i];
int orderId = chromOrder[valId];
if (orderId == -1)
continue;
int y = orderId * h;
safef(pixelStr, sizeof(pixelStr), "%d,%d", x, y);
struct hashEl *pixelEl = hashLookup(pixelHash, pixelStr);
if (!pixelEl)
{
hm = AllocA(struct hmPixel);
hm->x = x;
hm->y = y;
hm->w = w;
hm->h = h;
hm->val = 0.0;
hm->count = 0;
slAddHead(&hmList, hm);
hashAdd(pixelHash, pixelStr, hm);
}
else
hm = pixelEl->val;
if (isSNP)
{
if (val > 0.0)
{
hm->val += 1.0;
hm->count += 1;
}
}
else
{
hm->val += val;
hm->count += 1;
}
}
}
for (hm = hmList; hm ; hm = hm->next)
{
val = hm->val / (double) hm->count;
if(val > 0)
absVal = val;
else
absVal = -val;
absVal = absVal *gain;
int colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0)
colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if(val >= 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
vgBox(vg, hm->x, hm->y, hm->w, hm->h, valCol);
}
vgUnclip(vg);
freeHash(&pixelHash);
slFreeList(&hmList);
hmList = NULL;
pixelHash = hashNew(0);
}
freeHash(&pixelHash);
}
void drawTufteTukeyBox(struct vGfx *vg, struct heatmapLay *hl,
struct hmPixel *hm, double colorScale,
Color *upShades, Color *downShades, double absMax, int height)
/* Draws a dense, high data-to-ink ratio Tukey Boxplot, conceived by Edward Tufte
* Dot for the median, and then lines connecting the quartile to the whisker */
{
double val, absVal;
Color valCol;
/* Line connecting low quartile to low whisker, colored according to their average */
val = (hm->lowW + hm->lowQ)/2.0;
absVal = fabs(val);
int colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if(val > 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
val = hm->lowQ/absMax * height;
int lQ = height - round(val); // Low Quartile
val = hm->lowW/absMax * height;
int lW = height - round(val); // Low Whisker
int h = lW - lQ;
vgBox(vg, hm->x + hm->w/2, lQ, 1, h, valCol);
/* Line connecting upper quartile to upper whisker, colored according to their average */
val = (hm->highW + hm->highQ)/2.0;
absVal = fabs(val);
colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if(val > 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
val = hm->highQ/absMax * height;
int hQ = height - round(val); // High Quartile
val = hm->highW/absMax * height;
int hW = height - round(val); // High Whisker
h = hQ - hW;
vgBox(vg, hm->x + hm->w/2, hW, 1, h, valCol);
/* Draw a dot for the median */
val = hm->median/absMax * height;
int median = height - round(val);
vgBox(vg, hm->x + hm->w/2, median, 1, 1, MG_BLACK);
/* For scaling, check high and low whisker points */
if (hm->highW > hl->maxVal)
{
hl->maxVal = hm->highW;
hl->maxValPos = hW;
}
if (hm->lowW < hl->minVal)
{
hl->minVal = hm->lowW;
hl->minValPos = lW;
}
}
void drawStandardTukeyBox(struct vGfx *vg, struct heatmapLay *hl,
struct hmPixel *hm, double colorScale,
Color *upShades, Color *downShades, double absMax, int height)
/* Standard Tukey Boxplot, a box from upper->lower quartile, lines to whiskers, and a line
* for the median
*/
{
Color valCol;
double val = hm->median/absMax * height;
double absVal = fabs(hm->median);
/* Median a line, with 1 pixel buffer*/
int median = height - round(val);
vgBox(vg, hm->x+1, median, hm->w-2, 1, MG_BLACK); //valCol);
/* Lower whisker line */
val = (hm->lowW + hm->lowQ)/2.0;
absVal = fabs(val);
int colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if(val > 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
val = hm->lowQ/absMax * height;
int lQ = height - round(val); // Low Quartile
val = hm->lowW/absMax * height;
int lW = height - round(val); // Low Whisker
int h = lW - lQ;
vgBox(vg, hm->x + hm->w/2, lQ, 1, h, valCol);
vgBox(vg, hm->x+1, lW, hm->w-2, 1, valCol);
/* Upper whisker line*/
val = (hm->highW + hm->highQ)/2.0;
absVal = fabs(val);
colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if(val > 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
val = hm->highQ/absMax * height;
int hQ = height - round(val); // High Quartile
val = hm->highW/absMax * height;
int hW = height - round(val); // High Whisker
h = hQ - hW;
vgBox(vg, hm->x + hm->w/2, hW, 1, h, valCol);
vgBox(vg, hm->x+1, hW, hm->w-2, 1, valCol);
/* Upper/Lower Quartile box */
vgBox(vg, hm->x+1, lQ, hm->w-2, 1, MG_BLACK);
vgBox(vg, hm->x+1, hQ, hm->w-2, 1, MG_BLACK);
vgBox(vg, hm->x+1, hQ, 1, lQ - hQ, MG_BLACK);
vgBox(vg, hm->x + hm->w-2, hQ, 1, lQ - hQ, MG_BLACK);
/* For scaling, check high and low whisker points */
if (hm->highW > hl->maxVal)
{
hl->maxVal = hm->highW;
hl->maxValPos = hW;
}
if (hm->lowW < hl->minVal)
{
hl->minVal = hm->lowW;
hl->minValPos = lW;
}
}
double hmPixelComputeTukeySummary(struct genoHeatmap *gh, struct hmPixel *hmList)
/* Summary of data -- for drawing Tukey Box plot, returns absolute max for plotting purposes */
{
float gain = gh->gainFull;
struct hmPixel *hm;
double percentile = 0.25;
double absMax = INT_MIN;
for (hm = hmList; hm; hm = hm->next)
{
if (!hm->valList)
continue;
double count = (double) slCount(hm->valList);
slSort(&hm->valList, slDoubleCmp);
struct slDouble *sld = slElementFromIx(hm->valList, floor(percentile * count));
hm->lowQ = sld->val * gain;
sld = slElementFromIx(hm->valList, ceil((1.0-percentile) * (count - 1.0)));
hm->highQ = sld->val * gain;
hm->median = gain * slDoubleMedian(hm->valList);
double iqr = hm->highQ - hm->lowQ;
double lowBound = hm->lowQ - 1.5*iqr;
double highBound = hm->highQ + 1.5*iqr;
double tmp;
hm->lowW = INT_MAX;
hm->highW = INT_MIN;
for (sld = hm->valList; sld; sld = sld->next)
{
tmp = sld->val*gain;
if (tmp > lowBound && tmp < hm->lowW)
hm->lowW = tmp;
if (tmp < highBound && tmp > hm->highW)
hm->highW = tmp;
}
tmp = -1;
if (hm->lowW != INT_MAX && hm->highW != INT_MIN)
tmp = max(fabs(hm->lowW), fabs(hm->highW));
else if (hm->lowW != INT_MAX)
tmp = fabs(hm->lowW);
else if (hm->highW != INT_MIN)
tmp = fabs(hm->highW);
if (tmp > absMax)
absMax = tmp;
}
return absMax;
}
void drawTukeyBoxPlot(struct vGfx *vg, struct heatmapLay *hl, struct genoHeatmap *gh,
struct hmPixel *hmList, int width, int height)
{
double md = maxDeviation(gh->name);
double colorScale = COLOR_SCALE / md;
Color upShades[EXPR_DATA_SHADES];
Color downShades[EXPR_DATA_SHADES];
// default, expression data color scheme:
vgMakeColorGradient(vg, hl->zeroColor, hl->highColor, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, hl->zeroColor, hl->lowColor, EXPR_DATA_SHADES, downShades);
double absMax = hmPixelComputeTukeySummary(gh, hmList);
struct hmPixel *hm;
int midHeight = height/2;
vgSetClip(vg, 0, 0, width, height);
//vgBox(vg, 0, 0, width, height, MG_WHITE);
for (hm = hmList; hm ; hm = hm->next)
{
if (!hm->valList)
continue;
if (hm->w < 5) // dense version of the box, courtesy of Edward Tufte.
drawTufteTukeyBox(vg, hl, hm, colorScale, upShades, downShades, absMax, midHeight);
else // standard box plot
drawStandardTukeyBox(vg, hl, hm, colorScale, upShades, downShades, absMax, midHeight);
}
vgUnclip(vg);
}
void drawGenesetAverageLine(struct vGfx *vg, struct hmPixel *hmList,
double absMax, int height)
/* Draws a bargraph, specific for SNP data (for now) */
{
double absVal = fabs(absMax);
struct hmPixel *hm;
for (hm = hmList; hm; hm = hm->next)
{
int h = round(hm->avg/absVal * (double) height );
if (h == 0)
h = 1;
vgBox(vg, hm->x, height - h, hm->w, 2, MG_GRAY);
}
}
void drawBarGraph(struct vGfx *vg, struct heatmapLay *hl,
struct hmPixel *hm, double colorScale,
Color *upShades, Color *downShades, double absMax, int height)
/* Draws a bargraph, specific for SNP data (for now) */
{
Color valCol;
double absVal = fabs(absMax);
int colorIndex = EXPR_DATA_SHADES-1;
if(hm->avg > 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
int h = round(hm->avg/absVal * (double) height );
if (h == 0)
h = 1;
vgBox(vg, hm->x, height - h, hm->w, h, valCol);
/* For scaling, check high and low whisker points */
if (hm->avg > hl->maxVal)
{
hl->maxVal = hm->avg;
hl->maxValPos = height - h;
}
if (hm->avg < hl->minVal)
{
hl->minVal = 0.0; // Always set minimum val to zero for bar graphs.
hl->minValPos = height;
}
}
double hmPixelComputeBarGraphSummary(struct genoHeatmap *gh, struct hmPixel *hmList)
/* Summary of data -- Bar graph, currently specific for SNP bed15 files */
{
struct hmPixel *hm;
double absMax = INT_MIN;
for (hm = hmList; hm; hm = hm->next)
{
if (!hm->valList)
continue;
double sum = 0.0;
struct slDouble *sld;
for (sld = hm->valList; sld; sld = sld->next)
sum += sld->val;
hm->sum = sum;
hm->avg = sum / (double) hm->count;
if (hm->avg > absMax)
absMax = hm->avg;
}
if (absMax == 0.0)
absMax = 0.000000001;
return absMax;
}
void drawBarGraphPlot(struct vGfx *vg, struct heatmapLay *hl, struct genoHeatmap *gh,
struct hmPixel *hmList, int width, int height)
{ /* Specifically for SNP mutation data */
double md = maxDeviation(gh->name);
double colorScale = COLOR_SCALE / md;
Color upShades[EXPR_DATA_SHADES];
Color downShades[EXPR_DATA_SHADES];
// default, expression data color scheme:
vgMakeColorGradient(vg, hl->zeroColor, hl->highColor, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, hl->zeroColor, hl->lowColor, EXPR_DATA_SHADES, downShades);
double absMax = hmPixelComputeBarGraphSummary(gh, hmList);
struct hmPixel *hm;
vgSetClip(vg, 0, 0, width, height);
for (hm = hmList; hm ; hm = hm->next)
{
if (!hm->valList)
continue;
drawBarGraph(vg, hl, hm, colorScale, upShades, downShades, absMax, height);
}
vgUnclip(vg);
}
void drawGeneSetSNPSummary(struct vGfx *vg, struct genoHeatmap *gh,
struct heatmapLay *hlList, struct slName *sampleList, int height)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
if (!hlList)
return;
struct bed *nb = NULL;
double val;
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hmPixel *hm, *hmList = NULL;
struct hmPixel *hmGs, *hmGsList = NULL;
struct slName *sl;
struct slInt *si, *siList = NULL;
for (sl = sampleList; sl; sl = sl->next)
{
int i = slNameFindIx(gh->sampleList, sl->name);
struct slName *sample = slElementFromIx(gh->sampleList, i);
if (i < 0)
continue;
if (!slNameInList(sampleList, sample->name))
continue;
si = slIntNew(i);
slAddHead(&siList, si);
}
slReverse(&siList);
int sampleCount = slCount(siList);
int maxX = 0;
struct heatmapLay *hl;
for (hl = hlList; hl ; hl = hl->next)
{
struct hmElement *hEl;
hmGs = AllocA(struct hmPixel);
hmGs->x = hl->pixelStart;
hmGs->y = 0;
hmGs->w = hl->pixelEnd - hl->pixelStart;
hmGs->h = 0;
hmGs->val = 0.0;
hmGs->valList = NULL;
hmGs->count = 0;
hmGs->avg = -1.0;
slAddHead(&hmGsList, hmGs);
struct hash *txLength = hashNew(0);
for (hEl = hl->elements; hEl ; hEl = hEl->next)
{
nb = hEl->bed;
if (!nb)
continue;
int x = hEl->pixelStart;
int w = hEl->pixelEnd - hEl->pixelStart;
if (hEl->pixelEnd > maxX)
maxX = hEl->pixelEnd;
safef(pixelStr, sizeof(pixelStr), "%d", x);
struct hashEl *pixelEl = hashLookup(pixelHash, pixelStr);
if (!pixelEl)
{
hm = AllocA(struct hmPixel);
hm->x = x;
hm->y = 0;
hm->w = w;
hm->h = 0;
hm->val = 0.0;
hm->valList = NULL;
hm->count = 0;
slAddHead(&hmList, hm);
hashAdd(pixelHash, pixelStr, hm);
}
else
hm = pixelEl->val;
struct hashEl *el = hashLookup(txLength, hEl->name);
if (!el)
{
/* we're using score to store transcript length for SNP data */
int score = nb->score;
if (score == 0)
score = 1;
hashAddInt(txLength, hEl->name, score);
hmGs->count += score * sampleCount;
hm->count += score * sampleCount;
}
for(si = siList; si; si = si->next)
{
int i = si->val;
val = nb->expScores[i];
/* Count observed mutations */
hmGs->val += val;
struct slDouble *sld = slDoubleNew(val);
slAddHead(&hm->valList, sld);
}
}
freeHash(&txLength);
if (hmGs->count)
hmGs->avg = hmGs->val / (double) hmGs->count;
}
slReverse(&hmGsList);
drawBarGraphPlot(vg, hlList, gh, hmList, maxX, height);
drawGenesetAverageLine(vg, hmGsList, hlList->maxVal, height);
}
void drawGeneSetGenomicSummary(struct vGfx *vg, struct genoHeatmap *gh,
struct heatmapLay *hlList, struct slName *sampleList, int height)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
if (!hlList)
return;
struct bed *nb = NULL;
double val;
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hmPixel *hm, *hmList = NULL;
struct slName *sl;
struct slInt *si, *siList = NULL;
for (sl = sampleList; sl; sl = sl->next)
{
int i = slNameFindIx(gh->sampleList, sl->name);
struct slName *sample = slElementFromIx(gh->sampleList, i);
if (i < 0)
continue;
if (!slNameInList(sampleList, sample->name))
continue;
si = slIntNew(i);
slAddHead(&siList, si);
}
slReverse(&siList);
int maxX = 0;
struct heatmapLay *hl;
for (hl = hlList; hl ; hl = hl->next)
{
struct hmElement *hEl;
for (hEl = hl->elements; hEl ; hEl = hEl->next)
{
nb = hEl->bed;
if (!nb)
continue;
int x = hEl->pixelStart;
int w = hEl->pixelEnd - hEl->pixelStart;
if (hEl->pixelEnd > maxX)
maxX = hEl->pixelEnd;
safef(pixelStr, sizeof(pixelStr), "%d", x);
struct hashEl *pixelEl = hashLookup(pixelHash, pixelStr);
if (!pixelEl)
{
hm = AllocA(struct hmPixel);
hm->x = x;
hm->y = 0;
hm->w = w;
hm->h = 0;
hm->val = 0.0;
hm->valList = NULL;
hm->count = 0;
slAddHead(&hmList, hm);
hashAdd(pixelHash, pixelStr, hm);
}
else
hm = pixelEl->val;
for(si = siList; si; si = si->next)
{
int i = si->val;
val = nb->expScores[i];
struct slDouble *sld = slDoubleNew(val);
slAddHead(&hm->valList, sld);
}
}
}
drawTukeyBoxPlot(vg, hlList, gh, hmList, maxX, height);
}
void drawGeneSetSummary(struct vGfx *vg, char* database, char *chromHeatmap,
struct heatmapLay *hlList, struct slName *sampleList, int height)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
if (!hlList)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
if (sameString(gh->platform, "SNP"))
drawGeneSetSNPSummary(vg, gh, hlList, sampleList, height);
else
drawGeneSetGenomicSummary(vg, gh, hlList, sampleList, height);
}
void drawBedGraph(struct vGfx *vg, char* database,
struct heatmapLay *hl, char *tableName)
{
if (!hl)
return;
struct hashEl *el = hashLookup(ghHash, tableName);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", tableName);
struct bed *bg=NULL, *ibg=NULL;
float val;
double pixelsPerBase = 1.0/hl->basesPerPixel;
int x1,x2;
int y1,y2;
int w, h;
int start,end;
double gMin = chromGraphMin(tableName);
double gMax = chromGraphMax(tableName);
double height = hl->height;
double gScale = height/(gMax-gMin);
struct rgbColor color = chromGraphColor(tableName);
struct rgbColor black = {0,0,0};
Color valCol;
double md = 10.0;
double colorScale = COLOR_SCALE / md;
int nField = 4;
Color shadesOfColor[EXPR_DATA_SHADES];
if (sameWord(gh->dataType, "bed 4"))
vgMakeColorGradient(vg, &black, &color, 1, shadesOfColor);
else if (sameWord(gh->dataType, "bed 5"))
{
vgMakeColorGradient(vg, &black, &color, EXPR_DATA_SHADES, shadesOfColor);
nField = 5;
}
struct hmElement *hEl;
for(hEl = hl->elements; hEl; hEl = hEl->next)
{
bg = getBedGraph(gh, hEl->name, nField);
if (bg== NULL)
continue;
int chromX = hEl->pixelStart, chromY = 0;
int width = hEl->pixelEnd - hEl->pixelStart;
vgSetClip(vg, chromX, chromY, width, height);
/* Draw rest of points, connecting with line to previous point
* if not too far off. */
for(ibg = bg; ibg ; ibg = ibg->next)
{
start = ibg->chromStart;
end = ibg->chromEnd;
val = sqlFloat(ibg->name);
valCol = shadesOfColor[0];
if (sameWord(gh->dataType, "bed 5"))
{
int prob = ibg->score/md;
int colorIndex = abs(prob) * (EXPR_DATA_SHADES-1.0) * colorScale;
if (colorIndex < 0)
colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
valCol = shadesOfColor[colorIndex];
}
/* within a boundary */
if( val > gMax)
val = gMax;
else if (val < gMin)
val = gMin;
x1 = pixelsPerBase*start + chromX;
x2 = pixelsPerBase*end + chromX;
w = x2-x1;
if (val>0)
{
y1 = (height - (val - gMin) * gScale) + chromY; //smaller value
y2 = (height - (gMin -gMin) * gScale) + chromY;//greater value
h = y2-y1;
/* sometimes if the bed segment is too small, resulting in width=0,
* reset to 1 pixel length */
if ((w ==0) && (h != 0))
w=1;
vgBox(vg,x1,y1, w, h,valCol);
}
if (val<0)
{
y1 = (height - (val - gMax - gMin) *gScale) + chromY; //greater value
y2 = (height - (gMax - gMax - gMin)*gScale) + chromY; //smaller value
h = y1-y2;
/* sometimes if the bed segment is too small, resulting in width=0,
*reset to 1 pixel length */
if ((w ==0) && (h != 0))
w=1;
vgBox(vg,x1,y2,w,h,valCol);
}
}
vgUnclip(vg);
bedFree(&bg);
}
}
void drawChromSNPSummary(struct vGfx *vg, struct genoHeatmap *gh, struct heatmapLay *hl,
struct slName *sampleList, int height)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
struct bed *ghBed=NULL, *nb=NULL;
double pixelsPerBase = 1.0/hl->basesPerPixel;
double val;
int start, end;
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hash *txLength = hashNew(0);
struct hmPixel *hm, *hmList = NULL;
struct slName *sl;
struct slInt *si, *siList = NULL;
for (sl = sampleList; sl; sl = sl->next)
{
int i = slNameFindIx(gh->sampleList, sl->name);
if (i < 0)
continue;
si = slIntNew(i);
slAddHead(&siList, si);
}
slReverse(&siList);
int sampleCount = slCount(siList);
struct hmElement *hEl;
for(hEl = hl->elements; hEl; hEl = hEl->next)
{
ghBed = getChromHeatmapRange(gh, hEl->name, hEl->baseStart, hEl->baseEnd, FALSE);
int chromX = hEl->pixelStart;
for(nb = ghBed; nb; nb = nb->next)
{
start = nb->chromStart;
end = nb->chromEnd;
if (start < hEl->baseStart && end > hEl->baseEnd)
{ // Window is *inside* of a probe
start = hEl->baseStart;
end = hEl->baseEnd;
int width = hEl->pixelEnd - hEl->pixelStart;
pixelsPerBase = (double) width / (double) (end - start);
}
int x = pixelsPerBase * (start - hEl->baseStart) + chromX;
int w = pixelsPerBase * (end - start);
if (w == 0)
w = 1;
safef(pixelStr, sizeof(pixelStr), "%d", x);
struct hashEl *el = hashLookup(pixelHash, pixelStr);
if (!el)
{
hm = AllocA(struct hmPixel);
hm->x = x;
hm->y = 0;
hm->w = w;
hm->h = 0;
hm->val = 0.0;
hm->count = 0;
hm->valList = NULL;
slAddHead(&hmList, hm);
hashAdd(pixelHash, pixelStr, hm);
}
else
hm = el->val;
/* total kludge to get number of sites in a pixel in chromosome view
* assumes that a pixel doesn't have two identically sized genes */
int score = nb->score;
if (score == 0)
score = 1;
safef(pixelStr, sizeof(pixelStr), "%d,%d", x, score);
el = hashLookup(txLength, pixelStr);
if (!el)
{
hm->count += score * sampleCount;
hashAddInt(txLength, pixelStr, score);
}
for(si = siList; si; si = si->next)
{
int i = si->val;
val = nb->expScores[i];
struct slDouble *sld = slDoubleNew(val);
slAddHead(&hm->valList, sld);
}
}
bedFreeList(&ghBed);
}
drawBarGraphPlot(vg, hl, gh, hmList, hl->width, height);
}
void drawChromGenomicSummary(struct vGfx *vg, struct genoHeatmap *gh, struct heatmapLay *hl,
struct slName *sampleList, int height)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
struct bed *ghBed=NULL, *nb=NULL;
double pixelsPerBase = 1.0/hl->basesPerPixel;
double val;
int start, end;
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hmPixel *hm, *hmList = NULL;
struct slName *sl;
struct slInt *si, *siList = NULL;
for (sl = sampleList; sl; sl = sl->next)
{
int i = slNameFindIx(gh->sampleList, sl->name);
if (i < 0)
continue;
si = slIntNew(i);
slAddHead(&siList, si);
}
slReverse(&siList);
boolean useAccessTable = TRUE;
if (singleChromImage(hl))
useAccessTable = FALSE;
struct hmElement *hEl;
for(hEl = hl->elements; hEl; hEl = hEl->next)
{
ghBed = getChromHeatmapRange(gh, hEl->name, hEl->baseStart, hEl->baseEnd, useAccessTable);
int chromX = hEl->pixelStart;
for(nb = ghBed; nb; nb = nb->next)
{
start = nb->chromStart;
end = nb->chromEnd;
if (start < hEl->baseStart && end > hEl->baseEnd)
{ // Window is *inside* of a probe
start = hEl->baseStart;
end = hEl->baseEnd;
int width = hEl->pixelEnd - hEl->pixelStart;
pixelsPerBase = (double) width / (double) (end - start);
}
int x = pixelsPerBase * (start - hEl->baseStart) + chromX;
int w = pixelsPerBase * (end - start);
if (w == 0)
w = 1;
safef(pixelStr, sizeof(pixelStr), "%d", x);
struct hashEl *el = hashLookup(pixelHash, pixelStr);
if (!el)
{
hm = AllocA(struct hmPixel);
hm->x = x;
hm->y = 0;
hm->w = w;
hm->h = 0;
hm->val = 0.0;
hm->count = 0;
hm->valList = NULL;
slAddHead(&hmList, hm);
hashAdd(pixelHash, pixelStr, hm);
}
else
hm = el->val;
for(si = siList; si; si = si->next)
{
int i = si->val;
val = nb->expScores[i];
struct slDouble *sld = slDoubleNew(val);
slAddHead(&hm->valList, sld);
}
}
bedFreeList(&ghBed);
}
drawTukeyBoxPlot(vg, hl, gh, hmList, hl->width, height);
}
void drawChromSummary(struct vGfx *vg, char* database, struct heatmapLay *hl,
char *chromHeatmap, struct slName *sampleList, int height)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
if (!hl)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
if (sameString(gh->platform, "SNP"))
drawChromSNPSummary(vg, gh, hl, sampleList, height);
else
drawChromGenomicSummary(vg, gh, hl, sampleList, height);
}
void drawBackgroundLines(struct vGfx *vg, int startX, int startY, int width, int height)
{
int x, spacing = 20;
Color lightBlue = vgFindColorIx(vg, 220, 220, 255);
vgSetClip(vg, startX, startY, width, height);
for (x = startX; x < (startX + width); x += spacing)
vgBox(vg, x, startY, 1, height, lightBlue);
vgUnclip(vg);
}
void drawChromHeatmapsByPixel(struct vGfx *vg, char* database,
struct heatmapLay *hl, char *chromHeatmap)
/* Draw chromosome graph on all chromosomes in layout at given
* y offset and height. */
{
if (!hl)
return;
struct hashEl *el = hashLookup(ghHash, chromHeatmap);
struct genoHeatmap *gh =NULL;
if (el)
gh = el->val;
else
errAbort("No heatmap %s", chromHeatmap);
int *chromOrder = getBedOrder(gh);
struct bed *ghBed=NULL, *nb=NULL;
double pixelsPerBase = 1.0/hl->basesPerPixel;
double md = maxDeviation(gh->name);
double colorScale = COLOR_SCALE / md;
double val;
double absVal;
int valId;
Color valCol;
int start, end;
float gain = gh->gainFull;
Color upShades[EXPR_DATA_SHADES];
Color downShades[EXPR_DATA_SHADES];
vgMakeColorGradient(vg, hl->zeroColor, hl->highColor, EXPR_DATA_SHADES, upShades);
vgMakeColorGradient(vg, hl->zeroColor, hl->lowColor, EXPR_DATA_SHADES, downShades);
char pixelStr[128];
struct hash *pixelHash = hashNew(0);
struct hmPixel *hm, *hmList = NULL;
boolean isSNP = FALSE;
if (sameString(gh->platform, "SNP"))
isSNP = TRUE;
boolean useAccessTable = TRUE;
if (singleChromImage(hl))
useAccessTable = FALSE;
int sampleHeight = hl->sampleHeight;
int totalHeight = sampleHeight * heatmapHeight(gh);
int buffer = PROBE_MARK_HEIGHT + PROBE_MARK_BUFFER;
struct hmElement *hEl;
for(hEl = hl->elements; hEl; hEl = hEl->next)
{
ghBed = getChromHeatmapRange(gh, hEl->name, hEl->baseStart, hEl->baseEnd, useAccessTable);
int chromX = hEl->pixelStart, chromY = 0;
int width = hEl->pixelEnd - hEl->pixelStart;
vgSetClip(vg, chromX, chromY, width, totalHeight + buffer);
Color color = vgFindColorIx(vg, hl->zeroColor->r, hl->zeroColor->g, hl->zeroColor->b);
vgBox(vg, chromX, chromY, width, totalHeight, color);
for(nb = ghBed; nb; nb = nb->next)
{
start = nb->chromStart;
end = nb->chromEnd;
if (start < hEl->baseStart && end > hEl->baseEnd)
{ // Window is *inside* of a probe
start = hEl->baseStart;
end = hEl->baseEnd;
pixelsPerBase = (double) width / (double) (end - start);
}
int i;
for(i = 0; i < nb->expCount; ++i)
{
val = nb->expScores[i];
valId = nb->expIds[i];
int orderId = chromOrder[valId];
if (orderId == -1)
continue;
int y = chromY + orderId * sampleHeight;
int x = pixelsPerBase * (start - hEl->baseStart) + chromX;
int w = pixelsPerBase * (end - start);
if (w == 0)
w = 1;
safef(pixelStr, sizeof(pixelStr), "%d,%d", x, y);
struct hashEl *el = hashLookup(pixelHash, pixelStr);
if (!el)
{
hm = AllocA(struct hmPixel);
hm->x = x;
hm->y = y;
hm->w = w;
hm->h = sampleHeight;
hm->val = 0.0;
hm->count = 0;
slAddHead(&hmList, hm);
hashAdd(pixelHash, pixelStr, hm);
}
else
hm = el->val;
if (isSNP)
{
if (val > 0)
{
hm->val += 1.0;
hm->count += 1;
}
}
else
{
hm->val += val;
hm->count += 1;
}
}
}
for (hm = hmList; hm ; hm = hm->next)
{
val = hm->val / hm->count;
absVal = fabs(val) * gain;
int colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
/* Clip color index to fit inside of array, since we may have brightened it. */
if (colorIndex < 0) colorIndex = 0;
if (colorIndex >= EXPR_DATA_SHADES)
colorIndex = EXPR_DATA_SHADES-1;
if(val > 0)
valCol = upShades[colorIndex];
else
valCol = downShades[colorIndex];
/* Draw probe data */
vgBox(vg, hm->x, hm->y, hm->w, hm->h, valCol);
/* Draw probe mark at bottom */
vgBox(vg, hm->x, hl->height+PROBE_MARK_BUFFER,
hm->w, PROBE_MARK_HEIGHT, MG_BLACK);
}
freeHash(&pixelHash);
pixelHash = hashNew(0);
slFreeList(&hmList);
hmList = NULL;
bedFreeList(&ghBed);
vgUnclip(vg);
}
}
char *genomeGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
struct hvGfx *vg;
char *db = gh->database;
int totalW = hgHeatmapDefaultPixWidth;
hl->height = heatmapHeight(gh) * hl->sampleHeight;
/* Include room for probe marks */
int totalH = hl->height + PROBE_MARK_BUFFER + PROBE_MARK_HEIGHT;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char *strToHash = cartSettingsString(hgh2Prefix, "genomeGif");
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
if (sameWord(gh->dataType,"bed 15"))
drawChromHeatmapsByPixel(vg->vg, db, hl, gh->name);
else if (sameWord(gh->dataType, "bed 4") || sameWord(gh->dataType, "bed 5"))
drawBedGraph(vg->vg, db, hl, gh->name);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *genomeSummaryGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh, struct slName *sampleList,
int totalH, int subset, int subsetNum)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
struct hvGfx *vg;
char *db = gh->database;
int totalW = hgHeatmapDefaultPixWidth;
totalH = hgHeatmapDefaultSummaryHeight; //totalH;
hl->height = totalH;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char modifier[128];
safef(modifier, sizeof(modifier), "genomeSummaryGif,heatmap%dof%d", subset, subsetNum);
char *strToHash = cartSettingsString(hgh2Prefix, modifier);
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
if (sameWord(gh->dataType,"bed 15"))
drawChromSummary(vg->vg, db, hl, gh->name, sampleList, totalH);
else if (sameWord(gh->dataType, "bed 4") || sameWord(gh->dataType, "bed 5"))
drawBedGraph(vg->vg, db, hl, gh->name);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *genomeStatsGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh, float colorCutoff)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
if (!gh->anaResult)
return NULL;
struct hvGfx *vg;
int totalW = hgHeatmapDefaultPixWidth;
int totalH = hghBed5Height;
hl->height = totalH;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char *strToHash = cartSettingsString(hgh2Prefix, "genomeStatsGif");
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
// drawBackgroundLines(vg->vg, 0, 0, totalW, totalH);
if (sameWord(gh->dataType,"bed 15"))
drawChromStats(vg->vg, gh->anaResult->stats, hl,
0, gh->anaResult->max, colorCutoff);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *genesetStatsGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh, float colorCutoff)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
if (!gh->anaResultHash)
return NULL;
struct hvGfx *vg;
int totalW = hgHeatmapDefaultPixWidth;
int totalH = hghBed5Height;
hl->height = totalH;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char *strToHash = cartSettingsString(hgh2Prefix, "genesetStatsGif");
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
if (sameWord(gh->dataType,"bed 15"))
drawGeneSetStats(vg->vg, gh->anaResultHash->hash, hl,
0, 3, colorCutoff, gh->name);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *genesetGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
struct hvGfx *vg;
int totalW = hgHeatmapDefaultPixWidth;
hl->height = heatmapHeight(gh) * hl->sampleHeight;
int totalH = hl->height;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char *strToHash = cartSettingsString(hgh2Prefix, "genesetGif");
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
drawLayoutLines(vg->vg, hl, totalW);
if (sameWord(gh->dataType,"bed 15"))
drawGeneSetHeatmapsByPixel(vg->vg, gh->database, gh->name, hl); //, geneHash);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *genesetSummaryGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh, struct slName *sampleList,
int totalH, int subset, int subsetNum)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
struct hvGfx *vg;
int totalW = hgHeatmapDefaultPixWidth;
totalH = hgHeatmapDefaultSummaryHeight;
hl->height = totalH;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char modifier[128];
safef(modifier, sizeof(modifier), "genomeSummaryGif,heatmap%dof%d", subset, subsetNum);
char *strToHash = cartSettingsString(hgh2Prefix, modifier);
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
drawLayoutLines(vg->vg, hl, totalW);
if (sameWord(gh->dataType,"bed 15"))
drawGeneSetSummary(vg->vg, gh->database, gh->name, hl, sampleList, totalH);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *featureGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh, char *tableName)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
struct hvGfx *vg;
int totalW = hl->width;
int totalH = heatmapHeight(gh) * hl->sampleHeight;
hl->height = totalH;
if (ifDrawFeatureLabel(hl))
totalH += hghFeatureLabel;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char *strToHash = cartSettingsString(hgh2Prefix, "featureGif");
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
/* draw feature sorter */
if (sameWord(gh->dataType,"bed 15"))
drawFeatures(vg->vg, tableName, hl, 0, FALSE);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *featureSummaryGif(struct sqlConnection *conn, struct heatmapLay *hl,
struct genoHeatmap *gh, char *tableName, struct slName *sampleList,
int buffer, int subset, int subsetNum)
/* Create genome GIF file and HT that includes it. */
{
if (!gh || !hl)
return NULL;
struct hvGfx *vg;
int totalW = hl->width;
int totalH = hgHeatmapDefaultSummaryHeight;
hl->height = totalH;
if (totalW * totalH == 0)
return NULL;
boolean isBottomSubgroup = (subset + 1 == subsetNum);
if (ifDrawFeatureLabel(hl) && isBottomSubgroup)
totalH += hghFeatureLabel;
struct tempName md5Tn;
char modifier[128];
safef(modifier, sizeof(modifier), "featureSummaryGif,%dof%d", subset, subsetNum);
char *strToHash = cartSettingsString(hgh2Prefix, modifier);
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
/* draw feature sorter */
if (sameWord(gh->dataType,"bed 15"))
{
if (subsetNum > 1)
drawFeatureSubgroupBar(vg->vg, sampleList, hl->sampleHeight, isBottomSubgroup);
drawFeatureSummary(vg->vg, tableName, hl, 0, FALSE, sampleList, isBottomSubgroup);
}
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
char *subgroupGif(struct genoHeatmap *gh, char *tableName, int width, int sampleHeight)
/* Create genome GIF file and HT that includes it. */
{
if (!gh)
return NULL;
struct hvGfx *vg;
int totalW = width;
int totalH = heatmapHeight(gh) * sampleHeight;
if (totalW * totalH == 0)
return NULL;
struct tempName md5Tn;
char *strToHash = cartSettingsString(hgh2Prefix, "subgroupGif");
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
/* draw subgroups */
if (sameWord(gh->dataType,"bed 15"))
drawSubgroups(vg->vg, tableName, totalW, totalH, sampleHeight);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
void drawSimpleScale(struct vGfx *vg, struct heatmapLay *hlList, int width, int height)
{
if (!hlList)
return;
int minValPos = 0, maxValPos = 0;
double minVal = INT_MAX;
double maxVal = INT_MIN;
struct heatmapLay *hl;
for (hl = hlList; hl; hl = hl->next)
{
if (hl->minVal < minVal)
{
minVal = hl->minVal;
minValPos = hl->minValPos;
}
if (hl->maxVal > maxVal)
{
maxVal = hl->maxVal;
maxValPos = hl->maxValPos;
}
}
if (minVal == INT_MAX || maxVal == INT_MIN) // not scaled.
return;
int dashWidth = 5;
int buffer = 2;
vgSetClip(vg, 0, 0, width, height);
vgBox(vg, 0, 0, width, height, MG_WHITE);
/* draw line on right side */
vgBox(vg, width-1, maxValPos, 1, minValPos - maxValPos, MG_BLACK);
int yZero = height / 2;
if (minVal != maxVal)
{
double slope = ((double) (maxValPos - minValPos))/(maxVal - minVal);
double tmp = slope * (0.0 - maxVal) + maxValPos;
yZero = ceil(tmp) - 1;
/* Put dash at zero */
vgBox(vg, width - dashWidth, yZero, dashWidth, 1, MG_BLACK);
}
/* Put dash at max value */
vgBox(vg, width - dashWidth, maxValPos, dashWidth, 1, MG_BLACK);
/* Label max point */
int textHeight = maxValPos - hlList->fontHeight/2;
if (textHeight > yZero - hlList->fontHeight)
textHeight = yZero - hlList->fontHeight;
if (textHeight < 0)
textHeight = 0;
char maxValStr[128];
if (fabs(maxVal) > 1000.0)
safef(maxValStr, sizeof(maxValStr), "%1.1e", maxVal);
else if (fabs(maxVal) > 100.0)
safef(maxValStr, sizeof(maxValStr), "%3.1f", maxVal);
else if (fabs(maxVal) > 10.0)
safef(maxValStr, sizeof(maxValStr), "%2.2f", maxVal);
else if (fabs(maxVal) > 0 && fabs(maxVal) < 0.001)
safef(maxValStr, sizeof(maxValStr), "%1.1e", maxVal);
else if (fabs(maxVal) == 0)
safef(maxValStr, sizeof(maxValStr), "0");
else
safef(maxValStr, sizeof(maxValStr), "%1.3f", maxVal);
vgTextRight(vg, 0, textHeight, width - (dashWidth+buffer), hlList->fontHeight,
MG_BLACK, hlList->font, maxValStr);
/* Put dash at min value */
vgBox(vg, width-dashWidth, minValPos-1, dashWidth, 1, MG_BLACK);
/* Label min point */
textHeight = (minValPos-1) - hlList->fontHeight/2;
if (textHeight < yZero + hlList->fontHeight)
textHeight = yZero + hlList->fontHeight;
if (textHeight > height - hlList->fontHeight)
textHeight = height - hlList->fontHeight;
char minValStr[128];
if (fabs(minVal) > 1000.0)
safef(minValStr, sizeof(minValStr), "%1.1e", minVal);
else if (fabs(minVal) > 100.0)
safef(minValStr, sizeof(minValStr), "%3.1f", minVal);
else if (fabs(minVal) > 10.0)
safef(minValStr, sizeof(minValStr), "%2.2f", minVal);
else if (fabs(minVal) > 0 && fabs(minVal) < 0.001)
safef(minValStr, sizeof(minValStr), "%1.1e", minVal);
else if (fabs(minVal) == 0.0)
safef(minValStr, sizeof(minValStr), "0");
else
safef(minValStr, sizeof(minValStr), "%1.3f", minVal);
vgTextRight(vg, 0, textHeight, width - (dashWidth+buffer), hlList->fontHeight,
MG_BLACK, hlList->font, minValStr);
vgUnclip(vg);
}
char *simpleScaleGif(struct heatmapLay *hl, struct genoHeatmap *gh, int subset, int subsetNum)
/* Create genome GIF file and HT that includes it. */
{
if (!gh)
return NULL;
int totalW = hgDefaultScalePixWidth;
int totalH = hl->height;
if (totalW * totalH == 0)
return NULL;
if (hl->maxVal == INT_MIN || hl->minVal == INT_MAX) // scaling not set.
return NULL;
struct tempName md5Tn;
char modifier[128];
safef(modifier, sizeof(modifier), "simpleScaleGif,%dof%d", subset, subsetNum);
char *strToHash = cartSettingsString(hgh2Prefix, modifier);
trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
{
struct hvGfx *vg = hvGfxOpenGif(totalW, totalH, md5Tn.forCgi);
drawSimpleScale(vg->vg, hl, totalW, totalH);
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}