7994cb46f31e9eea5b916d7c6b2f632de3d6cb47 kate Wed Feb 10 13:12:57 2016 -0800 Support for multiple GTEx releases/versions (to add V6 track). refs #15645 diff --git src/hg/hgc/gtexClick.c src/hg/hgc/gtexClick.c index 337d186..6a8eba3 100644 --- src/hg/hgc/gtexClick.c +++ src/hg/hgc/gtexClick.c @@ -1,203 +1,209 @@ /* Details pages for GTEx tracks */ /* Copyright (C) 2015 The Regents of the University of California * See README in this or parent directory for licensing information. */ #include "common.h" #include "hash.h" #include "hdb.h" #include "hvGfx.h" #include "trashDir.h" #include "hgc.h" #include "gtexGeneBed.h" #include "gtexTissue.h" #include "gtexSampleData.h" #include "gtexUi.h" +#include "gtexInfo.h" struct tissueSampleVals /* RPKM expression values for multiple samples */ { struct tissueSampleVals *next; char *name; /* GTEx tissue name */ char *description; /* GTEx tissue description */ int color; /* GTEx tissue color */ int count; /* number of samples */ double *vals; /* RPKM values */ double min, max; /* minimum, maximum value */ double q1, median, q3; /* quartiles */ struct slDouble *valList; /* used to create val array */ }; /********************************************************/ /* R implementation. Invokes R script */ void drawGtexRBoxplot(struct gtexGeneBed *gtexGene, struct tissueSampleVals *tsvList, - boolean doLogTransform) + boolean doLogTransform, char *version) /* Draw a box-and-whiskers plot from GTEx sample data, using R boxplot */ { /* Create R data frame. This is a tab-sep file, one row per sample, * with columns for sample, tissue, rpkm */ struct tempName dfTn; trashDirFile(&dfTn, "hgc", "gtexGene", ".df.txt"); FILE *f = fopen(dfTn.forCgi, "w"); if (f == NULL) errAbort("can't create temp file %s", dfTn.forCgi); fprintf(f, "sample\ttissue\trpkm\n"); struct tissueSampleVals *tsv; int sampleId=1; int i; for (tsv = tsvList; tsv != NULL; tsv = tsv->next) { int count = tsv->count; for (i=0; iname, tsv->vals[i]); } fclose(f); // Plot to PNG file struct tempName pngTn; trashDirFile(&pngTn, "hgc", "gtexGene", ".png"); char cmd[256]; /* Exec R in quiet mode, without reading/saving environment or workspace */ -safef(cmd, sizeof(cmd), "Rscript --vanilla --slave hgcData/gtexBoxplot.R %s %s %s %s", +safef(cmd, sizeof(cmd), "Rscript --vanilla --slave hgcData/gtexBoxplot.R %s %s %s %s %s", gtexGene->name, dfTn.forCgi, pngTn.forHtml, - doLogTransform ? "log=TRUE" : "log=FALSE"); + doLogTransform ? "log=TRUE" : "log=FALSE", version); int ret = system(cmd); if (ret == 0) { printf("
\n", pngTn.forHtml); //printf("
\n", //pngTn.forHtml, imageWidth, imageHeight); //pngTn.forHtml, 900, 500); } } struct gtexGeneBed *getGtexGene(char *item, char *table) /* Retrieve gene info for this item from the main track table */ { struct gtexGeneBed *gtexGene = NULL; struct sqlConnection *conn = hAllocConn(database); char **row; char query[512]; struct sqlResult *sr; if (sqlTableExists(conn, table)) { sqlSafef(query, sizeof(query), "select * from %s where name = '%s'", table, item); sr = sqlGetResult(conn, query); row = sqlNextRow(sr); if (row != NULL) { gtexGene = gtexGeneBedLoad(row); } sqlFreeResult(&sr); } hFreeConn(&conn); return gtexGene; } struct tissueSampleVals *getTissueSampleVals(struct gtexGeneBed *gtexGene, boolean doLogTransform, - double *maxValRet) + char *version, double *maxValRet) /* Get sample data for the gene. Optionally log10 it. Return maximum value seen */ { // TODO: support version table name. Likely move to lib. struct hash *tsHash = hashNew(0); struct tissueSampleVals *tsv; struct hashEl *hel; struct slDouble *val; double maxVal = 0; struct gtexSampleData *sd = NULL; char query[256]; char **row; +char buf[256]; char *sampleDataTable = "gtexSampleData"; +safef(buf, sizeof(buf), "%s%s", sampleDataTable, version); struct sqlConnection *conn = hAllocConn("hgFixed"); -assert(sqlTableExists(conn, sampleDataTable)); +assert(sqlTableExists(conn, buf)); sqlSafef(query, sizeof(query), "select * from %s where geneId='%s'", - sampleDataTable, gtexGene->geneId); + buf, gtexGene->geneId); struct sqlResult *sr = sqlGetResult(conn, query); while ((row = sqlNextRow(sr)) != NULL) { sd = gtexSampleDataLoad(row); if ((hel = hashLookup(tsHash, sd->tissue)) == NULL) { AllocVar(tsv); hashAdd(tsHash, sd->tissue, tsv); } else tsv = (struct tissueSampleVals *)hel->val; maxVal = max(maxVal, sd->score); val = slDoubleNew(sd->score); slAddHead(&tsv->valList, val); } /* Fill in tissue descriptions, fill values array and calculate stats for plotting Then make a list, suitable for sorting by tissue or score NOTE: Most of this not needed for R implementation */ struct gtexTissue *tis = NULL, *tissues = gtexGetTissues(); struct tissueSampleVals *tsList = NULL; int i; if (doLogTransform) maxVal = log10(maxVal+1.0); for (tis = tissues; tis != NULL; tis = tis->next) { tsv = hashFindVal(tsHash, tis->name); if (tsv == NULL) { /* no non-zero values for this tissue/gene */ AllocVar(tsv); val = slDoubleNew(0.0); slAddHead(&tsv->valList, val); } tsv->name = tis->name; tsv->description = tis->description; tsv->color = tis->color; int count = tsv->count = slCount(tsv->valList); double *vals = AllocArray(tsv->vals, count); for (i=0; ivalList); if (doLogTransform) vals[i] = log10(val->val+1.0); else vals[i] = val->val; } doubleBoxWhiskerCalc(tsv->count, tsv->vals, &tsv->min, &tsv->q1, &tsv->median, &tsv->q3, &tsv->max); slAddHead(&tsList, tsv); } if (maxValRet != NULL) *maxValRet = maxVal; return tsList; } void doGtexGeneExpr(struct trackDb *tdb, char *item) /* Details of GTEx gene expression item */ { struct gtexGeneBed *gtexGene = getGtexGene(item, tdb->table); if (gtexGene == NULL) errAbort("Can't find gene %s in GTEx gene table %s\n", item, tdb->table); genericHeader(tdb, item); // TODO: link to UCSC gene printf("Gene: %s
", gtexGene->name); char query[256]; +uglyf("transcriptId=%s\n", gtexGene->transcriptId); sqlSafef(query, sizeof(query), "select kgXref.description from kgXref, knownToEnsembl where knownToEnsembl.value='%s' and knownToEnsembl.name=kgXref.kgID", gtexGene->transcriptId); struct sqlConnection *conn = hAllocConn(database); char *desc = sqlQuickString(conn, query); hFreeConn(&conn); if (desc != NULL) printf("Description: %s
\n", desc); printf("Ensembl ID: %s
\n", gtexGene->geneId); printf("View at GTEx portal
\n", gtexGene->geneId); puts("

"); boolean doLogTransform = cartUsualBooleanClosestToHome(cart, tdb, FALSE, GTEX_LOG_TRANSFORM, GTEX_LOG_TRANSFORM_DEFAULT); double maxVal = 0.0; -struct tissueSampleVals *tsvs = getTissueSampleVals(gtexGene, doLogTransform, &maxVal); - -drawGtexRBoxplot(gtexGene, tsvs, doLogTransform); +char *versionSuffix = gtexVersionSuffix(tdb->table); +struct tissueSampleVals *tsvs = getTissueSampleVals(gtexGene, doLogTransform, + versionSuffix, &maxVal); +char *version = gtexVersion(tdb->table); +drawGtexRBoxplot(gtexGene, tsvs, doLogTransform, version); printTrackHtml(tdb); }