2012/10/09

Making Color Ramps in Matlab

When visualizing an array of data in a heatmap, a good color map makes a world of difference.

Thanks to my work in 'omics (i.e. transcriptomics - microarrays and RNASeq) I've looked at a lot of heatmaps over the past couple of years, and generated quite a few to boot.  Back in my Matlab heavy grad school days, I was generally happy with the default 'jet' color scheme (which given it's double rainbow-eseque aesthetics would make some individuals on this planet overly emotional).  Suffice it to say, I was a bit wary of straying far from the available maps (the others I used semi-regularly were "bone", "gray", and "hot").

Today I needed to create a nice color ramp in a GUI tool I've developed in Matlab for a dataset that spanned [-Inf, Inf].  Ideally, it should have three color stops:
  • a "cool" color for extreme negative values
  • a neutral color for 0
  • a "hot" color for extreme positive values
The most "viewable" ramp of this sort (e.g. one that the color non-blinds and color blinds can equally enjoy) would be:
  • blue
  • black
  • yellow

If I were generating this ramp in R it would be quite trivial with the colorRampPalette() function:
bky.ramp = colorRampPalette(c('blue', 'black', 'yellow'))

The above line would create a function bky.ramp() that you could use to specify a ramping palette of arbitrary length for a heatmap() (or any other plotting function):
heatmap(X, col=bky.ramp(256))


Doing this in Matlab is similar, but a tad more obscure.  If you look at the help for the colormap() function it says:
A colormap is an m-by-3 matrix of real numbers between 0.0 and 1.0. Each row is an RGB vector that defines one color. The kth row of the colormap defines the kth color, where map(k,:) = [r(k) g(k) b(k)]) specifies the intensity of red, green, and blue.
colormap(map) sets the colormap to the matrix map. If any values in map are outside the interval [0 1], you receive the error Colormap must have values in [0,1].
 I know that the colors I need are:
  • blue = [0 0 1]
  • black = [0 0 0]
  • yellow = [1 1 0]
but how do I ramp between them?  Well for that you need interp1():
 interp1 1-D interpolation (table lookup)
    YI = interp1(X,Y,XI) interpolates to find YI, the values of the
    underlying function Y at the points in the array XI. X must be a
    vector of length N.
    If Y is a vector, then it must also have length N, and YI is the
    same size as XI.  If Y is an array of size [N,D1,D2,...,Dk], then
    the interpolation is performed for each D1-by-D2-by-...-Dk value
    in Y(i,:,:,...,:).
    If XI is a vector of length M, then YI has size [M,D1,D2,...,Dk].
    If XI is an array of size [M1,M2,...,Mj], then YI is of size
    [M1,M2,...,Mj,D1,D2,...,Dk].
In its simplest invocation, it does linear interpolation between supplied points in Y over points XI. How is this used to create a BKY color ramp with 256 levels?  Like so:
bkyramp = interp1([blue; black; yellow], linspace(1,3,256));

If you're the type that likes to encapsulate things in reusable functions (which I am), you end up with something like this: