ngsplot heatmap scale discrepancy

Question

[Sadia]

Hi, I'm using ngsplot to generate a heatmap for ChIP-seq coverage. However, I would like to replot the heatmap myself with the underlying data provided by ngsplot.

I'm a bit confused about the scaling though. The heatmap generated by ngsplot has a scale of (0, 0.9) whereas the max value of the hm data provided by the program goes to 2.5. Is there a subset of data (outliers) that ngsplot ignores during the generation of its plots?

[Friederike]

which command are you using and did you check the help section/documentation of that function? Most likely, there's some pruning of extreme values so that the dynamic range is not dominated by outliers with very high or very low numbers, but there may also be some z-score normalization or the like going on.

[Friederike]

To follow up on this, this is what the plotheat function seems to do under the hood:

# Function to plot heatmaps.
plotheat <- functionplot.name, width, height, pointsize, 
    enrichList, title2plot, flood.q=.01){

    # Setup basic parameters.
    ncolor <- 256
    enrich.palette <- colorRampPalette(c('white', 'red'))
    np <- length(enrichList)

    # Filter genes with zero sd. --> this will influence the number of genes that are shown!
    for(i in 1:np){
        g.sd <- apply(enrichList[[i]], 1, sd)
        enrichList[[i]] <- enrichList[[i]][g.sd > 0, ]
    }

  ##... (removed some code)

    # Flooding extreme values which are identified by quantiles. --> this is probably what's 
   # influencing the color scheme the most! Note that flood.q is a variable that can be defined when running plotheat()
    flood.pts <- quantile(c(mat, recursive=T), c(flood.q, 1-flood.q))
    mat[mat < flood.pts[1]] <- flood.pts[1]
    mat[mat > flood.pts[2]] <- flood.pts[2]

##... more code

Code from here. Disclaimer: I've never used that package so there may be additional intricacies.

[Sadia]

I'm looking at the code, but I'm not able to tell exactly parameter is used as the cutoff limit. Could you shed some light on this?

The function I'm using within the program:

https://github.com/shenlab-sinai/ngsplot/blob/develop/bin/ngs.plot.r

[Friederike]

what's the command you're typing?

cut-off parameter: flood.q

Answer 1

[Friederike]

Based on lines 347-358 of ngsplot.r, you should obtain a file "heatmap.RData", which, if loaded via load("heatmap.RData"), should also tell you which flood.frac value was used (default is 0.1 for flood.q and 0.02 for flood.frac as far as I can see)

# Heatmap R data.
if(galaxy==1) {
     heat.dat <- file.path(oname1, 'heatmap.RData')
} else {
     heat.dat <- file.path(oname, 'heatmap.RData')
}
save(reg.list, uniq.reg, ng.list, pts, enrichList, v.low.cutoff, go.algo, 
     ctg.tbl, bam.pair, xticks, flood.frac, hm.color, unit.width, rr, 
     go.list, color.scale, v.lib.size, font.size, go.paras, low.count,
     color.distr, 
     file=heat.dat)
cat("Done\n")

[Sadia]

The command I used is the basic one:

ngs.plot.r -G mm9 -R bed -C configuration-file.txt -O outputfile.txt

I loaded heatmap.RData and it does have the parameters. So thank you for that!

Do you know what the difference between flood.q or flood.frac is? My file only contains flood.frac anyways and it is 0.02. So does that translate to 98th percentile of values that are included in the plot?

[Friederike]

ngs.plot.r -G mm9 -R bed -C configuration-file.txt -O outputfile.txt

I would guess that you can tune that setting in configuration-file.txt

Do you know what the difference between flood.q or flood.frac is

No.

So does that translate to 98th percentile of values that are included in the plot?

This is what the quantile function does:

> quantile(rnorm(100), c(0.02, 1-0.02))
       2%       98% 
-2.075636  2.147820 

[Friederike]

I've moved the comment about heatmap.RData to an answer since it seemed to help the most. Feel free to challenge that or, otherwise, accept the answer by clicking on the thumbs-up symbol.

[Sadia]

Thanks for your help!