R/html.R

Defines functions .pubmedLink .ncbiGeneLink .geneTable .getGeneAnno .geneReport .sortGeneTable .getHTMLOfMarkedPathway .makeGGraphHTML .makeKGraphHTML .makeKPathHTML .makeVolcHTML .makeHmapHTML .viewGraph .viewPath .viewSet .makeView .createIndex eaBrowse

Documented in eaBrowse

# 
# Author: ludwig geistlinger
# Date: June 24th 2010
#
# description:  collection of functions that convert outputs of
#               enrichment methods to html, i.e. visualizes them
#               user-friendly
#
# Update:
#   27 May 2015: 
#       - major reworking of HTML output based on ReportingTools & hwriter 
#       - interactive graphics based on imageMap from geneplotter & biocGraph
#
###############################################################################


#' Exploration of enrichment analysis results
#' 
#' Functions to extract a flat gene set ranking from an enrichment analysis
#' result object and to detailedly explore it.
#' 
#' 
#' @aliases eaBrowse ea.browse gsRanking gs.ranking
#' @param res Enrichment analysis result list (as returned by the functions
#' \code{\link{sbea}} and \code{\link{nbea}}).
#' @param nr.show Number of gene sets to show.  As default all statistically
#' significant gene sets are displayed.
#' @param graph.view Optional.  Should a graph-based summary (reports and
#' visualizes consistency of regulations) be created for the result?  If
#' specified, it needs to be a gene regulatory network, i.e. either an absolute
#' file path to a tabular file or a character matrix with exactly *THREE* cols;
#' 1st col = IDs of regulating genes; 2nd col = corresponding regulated genes;
#' 3rd col = regulation effect; Use '+' and '-' for activation/inhibition.
#' @param html.only Logical.  Should the html file only be written (without
#' opening the browser to view the result page)? Defaults to FALSE.
#' @param out.dir Output directory. If \code{NULL}, defaults to a 
#' timestamp-generated subdirectory of \code{configEBrowser("OUTDIR.DEFAULT")}. 
#' @param report.name Name of the HTML report. If \code{NULL}, defaults to the
#' enrichment method used.
#' @param signif.only Logical.  Display only those gene sets in the ranking,
#' which satisfy the significance level? Defaults to TRUE.
#' @return gsRanking: \code{\linkS4class{DataFrame}} with gene sets ranked by
#' the corresponding p-value;
#' 
#' eaBrowse: none, opens the browser to explore results.
#'
#' If not instructed otherwise (via argument \code{out.dir}), 
#' the main HTML report and associated files are written to 
#' \code{configEBrowser("OUTDIR.DEFAULT")}. 
#' See \code{?configEBrowser} to change the location. 
#' If \code{html.only=FALSE}, the HTML report will automatically be opened in 
#' your default browser.
#' @author Ludwig Geistlinger <Ludwig.Geistlinger@@sph.cuny.edu>
#' @seealso \code{\link{sbea}}, \code{\link{nbea}},
#' \code{\link{combResults}}
#' @examples
#' 
#'     
#'     # real data
#'     exprs.file <- system.file("extdata/exprs.tab", package="EnrichmentBrowser")
#'     cdat.file <- system.file("extdata/colData.tab", package="EnrichmentBrowser")
#'     rdat.file <- system.file("extdata/rowData.tab", package="EnrichmentBrowser")
#'     probeSE <- readSE(exprs.file, cdat.file, rdat.file)
#'     geneSE <- probe2gene(probeSE) 
#'     geneSE <- deAna(geneSE)
#'     metadata(geneSE)$annotation <- "hsa"
#' 
#'     # artificial enrichment analysis results
#'     gs <- makeExampleData(what="gs", gnames=names(geneSE))
#'     ea.res <- makeExampleData(what="ea.res", method="ora", se=geneSE, gs=gs)
#' 
#'     # (5) result visualization and exploration
#'     gsRanking(ea.res)
#'  
#'     out.dir <- configEBrowser("OUTDIR.DEFAULT") 
#'     eaBrowse(ea.res, out.dir=out.dir, report.name="oraReport")
#' 
#' @export eaBrowse
eaBrowse <- function(res, nr.show=-1, 
    graph.view=NULL, html.only=FALSE, out.dir=NULL, report.name=NULL)
{
    isAvailable("ReportingTools", type="software")

    method <- ifelse(is(res$method, "character"), res$method, NA)
    se <- res$se
    alpha <- res$alpha
    gs <- res$gs
   
    # create out dir
    if(is.null(out.dir)) 
    {
        out.dir <- configEBrowser("OUTDIR.DEFAULT")
        stamp <- format(Sys.time(), "%a_%b%d_%Y_%H%M%S")
        out.dir <- file.path(out.dir, stamp)
    }
    else out.dir <- path.expand(out.dir)
    if(!file.exists(out.dir)) dir.create(out.dir, recursive=TRUE)
    if(is.null(report.name)) report.name <- method
    
    # how many gene sets to show in the output?
    if(nr.show < 1) nr.show <- res$nr.sigs
    stopifnot(nr.show > 0)
    if(nr.show > nrow(res$res.tbl)) nr.show <- nrow(res$res.tbl)

    # add description & nr.genes per gene set
    res <- res$res.tbl[seq_len(nr.show),]

    # expecting transition from gene set lists to collections in the near future
    # gsc <- res$gsc
    gsc <- .gsList2Collect(gs[res[,1]])
    res[,1] <- vapply(res[,1], 
        function(s) unlist(strsplit(s, "_"))[1],
        character(1))
    
    is.kegg <- is(collectionType(gsc[[1]]), "KEGGCollection")
    is.go <- is(collectionType(gsc[[1]]), "GOCollection")

    gs.title <- vapply(gsc, description, character(1))
    nr.genes <- vapply(gsc, function(g) length(geneIds(g)), integer(1))

    cnames <- c(colnames(res)[1], "TITLE")
    resn <- DataFrame(res[,1], gs.title)
    if(!("NR.GENES" %in% colnames(res)))
    {
        cnames <- c(cnames, "NR.GENES")
        resn <- DataFrame(resn, nr.genes)
    }
    cnames <- c(cnames, colnames(res)[2:ncol(res)])
    resn <- DataFrame(resn, res[,2:ncol(res)]) 
    colnames(resn) <- cnames
    res <- resn

    # make gene pages
    im <- incidence(gsc)
    org <- organism(gsc[[1]])
    if(is.na(org)) org <- metadata(se)$annotation
    if(!length(org)) stop("Organism annotation not found!\n", 
        "Organism under study must be annotated via metadata(se)$annotation")

    message("Creating gene report ...")
    se <- se[colnames(im),]
    rcols <- sapply(c("FC.COL", "ADJP.COL"), configEBrowser)
    fDat <- rowData(se)[,rcols]
    fDat <- as.data.frame(fDat)
    gt <- suppressMessages(.geneTable(im, org, fcs=fDat))
    gn.cols <- sapply(c("SYM.COL", "GN.COL"), configEBrowser)
    rowData(se)[,gn.cols] <- DataFrame(gt[,gn.cols]) 
    gt.reps <- sapply(gsc, function(s) .geneReport(s, gt, out.dir)) 
    
    # make gene set page
    # (1) link .geneReport
    link <- paste0(names(gsc), ".html")
    res[,"NR.GENES"] <- hwriter::hwrite(res[,"NR.GENES"], link=link, table=FALSE)
   
    # set view: volcano, heatmap 
    message("Creating set view ...")
    out.prefix <- file.path(out.dir, names(gsc))
    names(out.prefix) <- names(gsc)
    vcol <- sapply(gsc, function(s) 
        .viewSet(se[geneIds(s),], out.prefix[setName(s)]))
    vcol <- hwriter::hwriteImage(sub("sview.html", "volc.png", vcol),
        link=vcol, table = FALSE, height=50, width=50, target="_blank")
    res <- DataFrame(res, vcol)
    colnames(res)[ncol(res)] <- "SET.VIEW" 

    # path view: kegg maps
    if(is.kegg)
    { 
        message("Creating kegg view ...")
        isAvailable("pathview", type="software")
        vcol <- sapply(gsc, function(s) 
            .viewPath(setName(s), se[geneIds(s),], out.prefix[setName(s)]))
        vcol <- hwriter::hwriteImage(sub("kview.html", "kpath.png", vcol),
            link=vcol, table = FALSE, height=50, width=50, target="_blank")
        res <- DataFrame(res, vcol)
        colnames(res)[ncol(res)] <- "PATH.VIEW"
    }

    # graph view: ggeaGraphs
    if(!is.null(graph.view)) 
    {
        message("Creating graph view ...")
        vcol <- sapply(gsc, function(s) .viewGraph(se[geneIds(s),], .queryGRN(
            geneIds(s), graph.view, index=FALSE), alpha, out.prefix[setName(s)]))
        vcol <- hwriter::hwriteImage(sub("html$", "png", vcol),
            link=vcol, table = FALSE, height=50, width=50, target="_blank")
        res <- DataFrame(res, vcol)
        colnames(res)[ncol(res)] <- "GRAPH.VIEW"
    }
 
    # (2) link KEGG / GO 
    link <- NULL
    GS.COL <- configEBrowser("GS.COL")
    if(is.kegg) link <- sapply(gsc, function(s) 
        .getHTMLOfMarkedPathway(setName(s), geneIds(s)[fDat[geneIds(s),2] < alpha]))
    else if(is.go) link <- paste0(configEBrowser("GO.SHOW.URL"), res[,GS.COL])
    if(!is.null(link)) res[,GS.COL] <- 
        hwriter::hwrite(res[,GS.COL], link=link, table=FALSE)

    htmlRep <- ReportingTools::HTMLReport(shortName=report.name,
        title=paste(toupper(method), configEBrowser("RESULT.TITLE"), sep=" - "),
        basePath=dirname(out.dir), reportDirectory=basename(out.dir))
    res <- as.data.frame(res)
    ReportingTools::publish(res, htmlRep) 
    rep <- ReportingTools::finish(htmlRep)
    if(!html.only && interactive())
    { 
        message(paste("Your output files are in", out.dir))
        message(paste0("HTML report: ", report.name, ".html"))
        # this can be removed upon release of RStudio v1.2 (Summer 2018)
        if(Sys.getenv("RSTUDIO") == "1") rep <- URLencode(rep)
        browseURL(rep)
    }
}

.createIndex <- function(meth, comb, out.dir, report.name="index")
{
    indexPage <- ReportingTools::HTMLReport(shortName = report.name,
                        title = 'EnrichmentBrowser: Index of Result Files',
                        basePath = dirname(out.dir), 
                        reportDirectory = basename(out.dir))

    vcol <- "global_sview.html"
    de.plot <- hwriter::hwriteImage(sub("sview.html", "volc.png", vcol),
        link=vcol, table = FALSE, height=200, width=200, target="_blank")

    de.tag <- hwriter::hwrite(de.plot, border=0, br=TRUE)
    ReportingTools::publish(de.tag, indexPage)
    delink <- ReportingTools::Link("DE Measures for each Gene", "de.txt")
    ReportingTools::publish(delink, indexPage)
    
    if(comb)
    {
        res.tag <- hwriter::hwrite("Combined Results", heading=4)
        ReportingTools::publish(res.tag, indexPage)
        ttlink <- ReportingTools::Link("Top Table", "comb.html")
        ReportingTools::publish(ttlink, indexPage)
        frlink <- ReportingTools::Link("Full Ranking", "comb.txt")
        ReportingTools::publish(frlink, indexPage) 
    }

    for(m in meth)
    {
        res.tag <- hwriter::hwrite(paste(toupper(m), "Results"), heading=4)
        ReportingTools::publish(res.tag, indexPage)
        ttlink <- ReportingTools::Link("Top Table", paste0(m, ".html"))
        ReportingTools::publish(ttlink, indexPage)
        frlink <- ReportingTools::Link("Full Ranking", paste0(m, ".txt"))
        ReportingTools::publish(frlink, indexPage) 
    }
    index <- ReportingTools::finish(indexPage)
    if(interactive())
    {
        message(paste("HTML report:", paste0(report.name, ".html")))
        # this can be removed upon release of RStudio v1.2 (Summer 2018)
        if(Sys.getenv("RSTUDIO") == "1") index <- URLencode(index)
        browseURL(index)
    }
}

.makeView <- function(html1, html2, gene.html.pos=c("bottom", "topright"))
{
    gene.html.pos <- match.arg(gene.html.pos)
    s <- unlist(strsplit(html1, "_"))[1]
    
    head <- hwriter::hmakeTag("head", hwriter::hmakeTag("title", s))
    html1.tag <- hwriter::hmakeTag("frame", name="volc", src=html1)
    html2.tag <- hwriter::hmakeTag("frame", name="hmap", src=html2)
    html3.tag <- hwriter::hmakeTag("frame", name="gene", scrolling="auto")
    
    if(gene.html.pos == "topright")
    {
        bkp <- html3.tag
        html3.tag <- html2.tag
        html2.tag <- bkp
    }

    f1.tag <- hwriter::hmakeTag("frameset", 
        paste0(sub("</frame>", "", c(html1.tag, html2.tag)), collapse=""), 
        cols=paste0(configEBrowser("PLOT.WIDTH") + 30, ",*"), border=0)
    html3.tag <- sub("</frame>", "", html3.tag)
    f2.tag <- hwriter::hmakeTag("frameset", paste(f1.tag, html3.tag),
        rows=paste0(configEBrowser("PLOT.HEIGHT") + 30, ",*"), border=0)
    cont <- hwriter::hmakeTag("html", paste0(head,f2.tag))
    return(cont)
}

.viewSet <- function(se, out.prefix)
{
    out.files <- paste(out.prefix, 
        c("sview.html", "volc.png", "hmap.png"), sep="_" )
    
    if(!file.exists(out.files[1]))
    {
        ##
        # 1 PLOT: heatmap & volcano
        ##
        volc.html <- .makeVolcHTML(se, out.files[2])
        hmap.html <- .makeHmapHTML(se, out.files[3])
        cont <- .makeView(volc.html, hmap.html) 
        cat(cont, file=out.files[1])
    }
    views <- basename(out.files[1])
    return(views)

    # flat file set report
    # Do we need this anymore?
#    report.file <- sub("pdf$", "txt", plot.file)
#    if(!file.exists(report.file))
#    {
#        fDat <- signif(fDat[order(fDat[,2]),], digits=2)
#        fDat <- cbind(rownames(fDat), fDat)
#        colnames(fDat)[1] <- "GENE"
#        write.table(fDat, file=report.file, row.names=FALSE, quote=FALSE, sep="\t")
#    }
#    views <- c(views, basename(report.file))
#    return(rev(views))
}

.viewPath <- function(s, se, out.prefix)
{
    org <- substring(s, 1, 3)
    pwy.id <- sub("^[a-z]{3}", "", s)
    fc <- rowData(se)[,configEBrowser("FC.COL")]
    gn.cols <- sapply(c("SYM.COL", "GN.COL"), configEBrowser)
    gnam <- apply(rowData(se)[,gn.cols], 1, paste, collapse=": ")
    names(fc) <- names(gnam) <- names(se)

    out.files <- paste(out.prefix, 
        c("kview.html", "kpath.png", "kgraph.png"), sep="_")
    
    OVIEW.IDS <- c("01100", "01110", "01120", "01130")
    if(!file.exists(out.files[1]) && !(pwy.id %in% OVIEW.IDS))
    {
        ##
        # 1 PLOT: kegg.native & kegg.graph
        ##
        kpath.html <- try(
            .makeKPathHTML(fc, pwy.id, org, out.files[2]),
            silent=TRUE
        )

        kgraph.html <- try(
            .makeKGraphHTML(fc, gnam, pwy.id, org, out.files[3]),
            silent=TRUE
        )
        cont <- .makeView(kgraph.html, kpath.html, gene.html.pos="topright") 
        cat(cont, file=out.files[1])
    }
    
    views <- basename(out.files[1])
    return(views)
} 

.viewGraph <- function(se, sgrn, alpha, out.prefix)
{
    out.files <- paste0(out.prefix, "_gview.", c("html", "png", "txt"))
    
    if(!file.exists(out.files[1]))
    {
        ##
        # 1 PLOT: ggeaGraph
        ##
        ggraph.html <- .makeGGraphHTML(se, sgrn, alpha, out.files[2])
        void.html <- "void.html"
        void.file <- file.path(dirname(out.prefix), void.html)
        cat(hwriter::hmakeTag('html'), file=void.file)
        cont <- .makeView(ggraph.html, void.html)   
        cat(cont, file=out.files[1])
    }
    
    views <- basename(out.files[1])
    return(views)
}

.makeHmapHTML <- function(se, img.file)
{
    width <- configEBrowser("PLOT.WIDTH") 
    height <- configEBrowser("PLOT.HEIGHT")

    # 1: make the plot
    # (a) heatmap 1: all features vs all samples
    expr <- assay(se)
    rownames(expr) <- rowData(se)[,configEBrowser("SYM.COL")]
    png(img.file, width=width, height=height)
    exprsHeatmap(expr=expr, grp=se[[configEBrowser("GRP.COL")]])
    dev.off()
    img.tag <- hwriter::hwriteImage(basename(img.file))

    # (b) heatmap 2: most signif features
    max.row <- 40
    fc <- abs(rowData(se)[,configEBrowser("FC.COL")])
    p <- -log(rowData(se)[,configEBrowser("ADJP.COL")], base=10)
    ind <- (fc >= 1) | (p >= 1)
    se <- se[ind,]
    if(nrow(se) > 1)
    {
        if(nrow(se) > max.row)
        {
            fc <- fc[ind]
            p <- p[ind]
            score <- sqrt(fc^2 + p^2)
            se <- se[order(score, decreasing=TRUE)[seq_len(max.row)],]
            #    # select most variable samples
            #if(ncol(se) > max.col)
            #{
            #    svar <- apply(assay(se), 2, var)
            #    se <- se[,order(svar, decreasing=TRUE)[seq_len(max.col)]]
            #}
        }
        expr <- assay(se)
        rownames(expr) <- rowData(se)[,configEBrowser("SYM.COL")]
        img2 <- sub(".png$", "2.png", img.file)
        png(img2, width=width, height=height)
        exprsHeatmap(expr=expr, grp=se[[configEBrowser("GRP.COL")]])
        dev.off()
        img.tag <- paste0(img.tag, hwriter::hwriteImage(basename(img2)))
    }
    
    # 2: make the html
    hmap.html <- sub("png$", "html", img.file)
    cont <- hwriter::hmakeTag('html', hwriter::hmakeTag('body', img.tag))
    cat(cont, file=hmap.html)
    return(basename(hmap.html))
}

.makeVolcHTML <- function(se, img.file)
{   
    width <- configEBrowser("PLOT.WIDTH") 
    height <- configEBrowser("PLOT.HEIGHT")

    fc <- rowData(se)[,configEBrowser("FC.COL")]
    p <- rowData(se)[,configEBrowser("ADJP.COL")]
    # 1: make the plot
    png(img.file, width=width, height=height)
        volcano(fc, p)
        mai.px <- par('mai') * 75
        usr.grd <- par('usr')
    dev.off()
    p <- -log(p, base=10)
    # 2: make the html

    # INFER LINK COORDS
	#
	# mai: a numerical vector of the form ‘c(bottom, left, top, right)
	#       which gives the margin size specified in inches
	#
	# default plot: (75px == 1 inch)
	#                       btm     left    top     right
	#       par('mai') :    1.02    0.82    0.82    0.42
	#           px:         76.5    61.5    61.5    31.5
	
	# actual plot region:
	#   - origin(x0,y0) on the top left is at (61.5, 61.5)
	#   - width: 500 - 61.5 - 31.5 = 407
	#   - height: 500 - 76.5 - 61.5 = 362
	
	plot.width <- width - sum(mai.px[c(2,4)])
	plot.height <- height - sum(mai.px[c(1,3)])
	
	# re-scale user grid
	min.x <- min(0, usr.grd[1])
	fc <- fc + abs(min.x)
	max.x <- usr.grd[2] + abs(min.x)
	min.y <- min(0, usr.grd[3])
	p <- p + abs(min.y)
	max.y <- usr.grd[4] + abs(min.y)
	
	y.scalef <- plot.height / max.y
	x.scalef <- plot.width / max.x
	
	cx <- fc * x.scalef + mai.px[2]
	cy <- p * y.scalef + mai.px[1]
	
	# turn y-coords around, (0,0) is topleft and not btmleft
	cy <- height - cy
	coord <- cbind(cx, cy, cx+10, cy-10)
        
    # volcano html
    volc.html <- sub("png$", "html", img.file) 
    con <- file(volc.html, open="w")
    
    gn.cols <- sapply(c("SYM.COL", "GN.COL"), configEBrowser)
    titles <- apply(rowData(se)[,gn.cols], 1, paste, collapse=": ") 
    refs <- paste0(configEBrowser("GENE.URL"), names(se))
    isAvailable("geneplotter", type = "software")
    geneplotter::imageMap(coord, con, 
        list(HREF=refs, TITLE=titles, TARGET=rep("gene", nrow(coord))), 
        basename(img.file))    
    close(con)
    return(basename(volc.html))
}

.makeKPathHTML <- function(fc, pwy.id, org, img.file)
{
    # 1: make the plot
    # run pathview for getting native overplotted image
    out.dir <- dirname(img.file)
    suppressWarnings(
        pathview::pathview(gene.data=fc, 
            pathway.id=pwy.id, species=org, kegg.dir=out.dir, out.suffix="kpath")
    )
    pv.out <- file.path(getwd(), paste0(org, pwy.id, ".kpath.png"))
    file.copy(from = pv.out, to = img.file)
    file.remove(pv.out)
    
    # 2: make the html
    kpath.html <- sub("png$", "html", img.file)
    img.tag <- hwriter::hwriteImage(basename(img.file))
    cont <- hwriter::hmakeTag('html', hwriter::hmakeTag('body', img.tag))
    cat(cont, file=kpath.html)
    return(basename(kpath.html))
}

.makeKGraphHTML <- function(fc, gname, pwy.id, org, img.file)
{
    width <- configEBrowser("PLOT.WIDTH") 
    height <- configEBrowser("PLOT.HEIGHT")

    # 1: make the plot
    # run pathview2 for getting the kegg graph
    out.dir <- dirname(img.file)
    png(img.file, width=width, height=height)
    par(mai=rep(0,4))
    suppressWarnings(suppressMessages(
        gr <- pathview2(gene.data=fc, 
            pathway.id=pwy.id, species=org, kegg.dir=out.dir)
    ))
    dev.off()   
 
    # 2: make the html
    kgraph.html <- sub("png$", "html", img.file)
    if(is(gr, "graph"))
    {
        nd <- nodeRenderInfo(gr)$kegg.ids
        nam <- sapply(names(nd), function(n) 
            ifelse(nd[[n]][1] %in% names(gname), 
                    gname[nd[[n]][1]], 
                    nodeRenderInfo(gr)$label[[n]]))
        names(nam) <- names(nd)
        kstr <- sapply(nd, function(n) 
            paste(paste(org, n, sep=":"), collapse="+"), USE.NAMES=FALSE)
        con <- file(kgraph.html, open="w")
        refs <- paste0(configEBrowser("KEGG.GENE.URL"), kstr)
        isAvailable("biocGraph", type = "software")
        biocGraph::imageMap(gr, con=con,
            tags=list(HREF=refs, TITLE = nam, TARGET = rep("gene", length(nd))),
            imgname=basename(img.file), width=width, height=height)    
        close(con)
    }
    else cat(hwriter::hmakeTag('html'), file=kgraph.html)
    return(basename(kgraph.html))
}

.makeGGraphHTML <- function(se, sgrn, alpha, img.file)
{

    width <- configEBrowser("PLOT.WIDTH") 
    height <- configEBrowser("PLOT.HEIGHT")

    sggeaGraph <- NULL
    if(nrow(sgrn) > 0)
        sggeaGraph <- .constructGGEAGraph(grn=sgrn, se=se, alpha=alpha)
    
    # txt report
    report.file <- sub("png$", "txt", img.file)
    if(!is.null(sggeaGraph))
    {
            consistency <- sggeaGraph@renderInfo@edges$label
            cons.tbl <- cbind(names(consistency), consistency)
            cons.tbl <- cons.tbl[order(as.numeric(consistency), decreasing=TRUE),]
            colnames(cons.tbl) <- c("EDGE", "CONSISTENCY")
            write.table(cons.tbl,
                file=report.file, row.names=FALSE, quote=FALSE, sep="\t")
    }
    else cat("No edges in network for this set!", file=report.file)

    # ggeaGraph png
    png(img.file, width=width, height=height)
    par(mai=rep(0,4))
    if(!is.null(sggeaGraph))
        sggeaGraph <- .plotGGEAGraph(sggeaGraph,
            show.scores=(numEdges(sggeaGraph) < configEBrowser("NR.SHOW")))
    else plot(NA, axes=FALSE, xlim=c(0,1), ylim=c(0,1),
        ylab="", xlab="", main="No edges in network for this set!")
    dev.off()
    
    ggraph.html <- sub("view.png$", "graph.html", img.file)
    if(!is.null(sggeaGraph))
    {
        gn.cols <- sapply(c("SYM.COL", "GN.COL"), configEBrowser)
        gnam <- apply(rowData(se)[,gn.cols], 1, paste, collapse=": ")
        names(gnam) <- names(se)

        # image map
        nd <- nodes(sggeaGraph)
        nam <- gnam[nd]
        con <- file(ggraph.html, open="w")
        refs <- paste0(configEBrowser("GENE.URL"),  nd)
        isAvailable("biocGraph", type = "software")
        biocGraph::imageMap(sggeaGraph, con=con,
            tags=list(HREF=refs, TITLE = nam, TARGET = rep("gene", length(nd))),
            imgname=basename(img.file), width=width, height=height)    
        close(con)
    }
    else cat(hwriter::hmakeTag('html'), file=ggraph.html)
    return(basename(ggraph.html))
}

.getHTMLOfMarkedPathway <- function(pwy, oids)
{
    pwy <- sub("^path:", "", pwy)
    oids <- gsub("[a-z]{3}:", "", oids)
    coids <- paste(oids, collapse="+")
    request <- pwy
    if(nchar(coids)) request <- paste(request, coids, sep="+")
    return(paste0(configEBrowser("KEGG.SHOW.URL"), request))
}

.sortGeneTable <- function(gt)
{
    gt <- as.data.frame(gt)
    num.cols <- sapply(gt, is.numeric)
    gt[,num.cols] <- signif(gt[,num.cols], digits=3)
    scols <- sapply(c("ADJP.COL", "FC.COL"), configEBrowser)
    sorting.df <- gt[, scols]
    sorting.df[,2] <- -abs(sorting.df[,2])
    gt <- gt[do.call(order, sorting.df), , drop=FALSE]
    return(gt)    
}

.geneReport <- function(s, gt, out.dir)
{
    # (0) extract gene information
    gt <- gt[geneIds(s),]
    gt <- .sortGeneTable(gt)
        
    # (1) html table
    htmlRep <- ReportingTools::HTMLReport(basePath=dirname(out.dir), 
        reportDirectory=basename(out.dir),
        shortName=setName(s), title=paste(setName(s), "Gene Report", sep=": "))
    ReportingTools::publish(gt, htmlRep, .modifyDF=list(.ncbiGeneLink))

    # (2) flat file
    fname <- paste0(setName(s), "_genes.txt")
    ofname <- file.path(out.dir, fname)
    if(!file.exists(ofname))
        write.table(gt, file=ofname, sep="\t", quote=FALSE, row.names=FALSE)
    dlink <- ReportingTools::Link("Download .txt", fname)
    ReportingTools::publish(dlink, htmlRep)

    rep <- ReportingTools::finish(htmlRep)
    return(rep)
}

.getGeneAnno <- function(ids, org, biotype=FALSE)
{
    # load org pkg
    org.pkg <- .org2pkg(org)
    isAvailable(org.pkg)
    org.pkg <- get(org.pkg)

    # (1) gene identification 
    EZ.COL <- configEBrowser("EZ.COL")
    gn.cols <- sapply(c("SYM.COL", "GN.COL"), configEBrowser)
    gt <- suppressMessages(AnnotationDbi::select(org.pkg, keys=ids,
            columns=gn.cols, keytype=EZ.COL))
	
	if(biotype)
	{
		useMart <- listDatasets <- useDataset <- getBM <- NULL
		isAvailable("biomaRt", type="software")
		
		id.type <- "entrezgene"
        message("Connecting to BioMart ...")
        ensembl <- useMart("ENSEMBL_MART_ENSEMBL", host = "www.ensembl.org")
        ds <- listDatasets(ensembl)[, "dataset"]
        ds <- grep(paste0("^", org), ds, value = TRUE)
        if (length(ds) == 0) 
            stop(paste("Mart not found for:", org))
        else if (length(ds) > 1) {
            message("Found several marts")
            sapply(ds, function(d) message(paste(which(ds == 
                d), d, sep = ": ")))
            n <- readline(paste0("Choose mart (1-", length(ds), 
                ") : "))
            ds <- ds[as.integer(n)]
        }
        ensembl <- useDataset(ds, mart = ensembl)
        attrs <- c(id.type, "gene_biotype")
        biot <- getBM(filters = id.type, attributes = attrs, 
            values = ids, mart = ensembl)
		ind <- match(ids, biot[,1])
		biot <- biot[ind,2]
		gt <- cbind(gt, biot)
		colnames(gt)[ncol(gt)] <- "BIOTYPE"
	}
    return(gt)
}

.geneTable <- function(im, org, fcs=NULL, grn=NULL)#, context="")
{
    # load org pkg
    org.pkg <- .org2pkg(org)
    isAvailable(org.pkg)
    org.pkg <- get(org.pkg)

    # (1) gene identification 
    EZ.COL <- configEBrowser("EZ.COL")
    gn.cols <- sapply(c("SYM.COL", "GN.COL"), configEBrowser)
    gt <- AnnotationDbi::select(org.pkg, 
            keys=colnames(im), columns=gn.cols, keytype=EZ.COL)

    # (2) fold change
    if(!is.null(fcs))
    {
        fcs[,1] <- round(fcs[,1], digits=2)
        fcs[,2] <- signif(fcs[,2], digits=2)
        gt <- cbind(gt, fcs) 
    }

    # (3) interactions
    if(!is.null(grn))
    {
        ias.per.gene <- sapply(colnames(im), 
            function(gene) grn[grn[,1] == gene | grn[,2] == gene,,drop=FALSE])
        nr.ias.per.gene <- sapply(ias.per.gene, nrow)
        gt <- cbind(gt, nr.ias.per.gene)
        colnames(gt)[ncol(gt)] <- "INTACTS"
    }

    # (4) nr.sets
#    gene.occ.freq <- colSums(im)
#    gt <- cbind(gt, gene.occ.freq)
#    colnames(gt)[ncol(gt)] <- "SETS"
#
#    # (5) pubmed
#    pmids <- AnnotationDbi::mapIds(org.pkg, keys=colnames(im), 
#        column="PMID", keytype=EZ.COL, multiVals="list")
#
#    # context search ?
#    # if(context != "") pmids <- search.abstracts(pmids, context=context)
#  
#    nr.articles <- sapply(pmids, length)
#    pubmed.urls <- sapply(pmids, function(p) paste0(
#        configEBrowser("PUBMED.URL"), paste(p, collapse=",")), USE.NAMES=FALSE)
#    articles <- paste(nr.articles, pubmed.urls)
#
#    gt <- cbind(gt, articles)
#    colnames(gt)[ncol(gt)] <- configEBrowser("PMID.COL")
    return(gt)
}


.ncbiGeneLink <- function(object, ...)
{
    EZ.COL <- configEBrowser("EZ.COL")
    col <- as.character(object[,EZ.COL])
    link <- paste0(configEBrowser("GENE.URL"), col)
    object[,EZ.COL] <- hwriter::hwrite(col, link=link, table=FALSE)
    return(object)
}

.pubmedLink <- function(object, ...)
{
    PMID.COL <- configEBrowser("PMID.COL")
    spl <- sapply(as.character(object[,PMID.COL]), 
        function(s) unlist(strsplit(s, " ")))
    spl <- t(spl)
    spl1 <- as.integer(spl[,1])
    object[,PMID.COL] <- hwriter::hwrite(spl1, link=spl[,2], table=FALSE)
    return(object)    
}

Try the EnrichmentBrowser package in your browser

Any scripts or data that you put into this service are public.

EnrichmentBrowser documentation built on Dec. 12, 2020, 2 a.m.