[R] densCols: what are the computed densities and how to create a legend

[Florian Hahne]

Hi Kate,
the function does not directly return the calculated density values, however those can be recomputed through the internal function .smoothScatterCalcDensity
This function is not exported from the grDevices name space, so you need to use the full qualifier:

dat <- cbind(rnorm(100), rnorm(100))
dc <- densCols(dat)
dd <- grDevices:::.smoothScatterCalcDensity(dat, nbin=128)

Densities are computed using a binning approach, and the grid of densities is available as the list item fhat. In order to map those to colors we need to find the absolute range of density values:
dens <- as.numeric(dd$fhat)
dens <- dens[dens>0]
range(dens)

And now we can construct a mapping of values to colors:
colLegend <- data.frame(density=seq(min(dens), max(dens), len=10), color=I(colorRampPalette(blues9[-(1:3)])(10)))

For plotting of a legend, you could try something like this:
plot(NA, xlim=c(0,10), ylim=c(0,11), type="n", ann=FALSE, axes=FALSE)
rect(0, 1:10, 1, 2:11, border=NA, col=colLegend$col)
text(2, (1:10)+0.5, signif(colLegend$density, 2))

And to get everything in one plot, something like this:
layout(matrix(rep(1:2, c(5,1)), nrow=1))
plot(dat, col=dc, pch=20)
par(mar=c(5,0,5,0))
plot(NA, xlim=c(0,10), ylim=c(0,11), type="n", ann=FALSE, axes=FALSE)
rect(0, 1:10, 1, 2:11, border=NA, col=colLegend$col)
text(2, (1:10)+0.5, signif(colLegend$density, 2), adj=0)

Please note that there might be some slight rounding errors here, since the colors for densCols are computed on a subset of bins (those in which the data points fall).  The bandwidth smoothing results in non-zero density values for neighboring bins, and since we map our colors to the overall range of densities, we can end up being slightly off. If you want this to be 100% accurate, you have to figure out the densities for the "filled" bins only and map color space to that, in which case you essentially reproduce all that densCols is doing (type densCols in your R session to see the source code and copy paste what you need...)

Hope that helps. 
Florian


On Apr 20, 2010, at 4:08 PM, Kate Zinszer wrote:

> 
> Hi Dr. Hahne,
> I saw that you were affiliated with this function so perhaps you have some insight in my problem that is stated below?
> 
> Many thanks,
> Kate
> 
> 
> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On Behalf Of Kate Zinszer
> Sent: 19 avril 2010 15:21
> To: r-help at R-project.org
> Subject: [R] densCols: what are the computed densities and how to create a legend
> 
> Hi,
> 
> I'm using the densCols function for a scatterplot and cannot figure out 1) how to extract the computed densities, and 2) how to create a legend based that represents the upper and lower ranges of the densities.
> 
> For example:
> 
> movers.den  <- densCols(move$x, move$y)
> table(movers.den)
> #08306B #083775 #083B7C #083D7E #3989C1 #3F8FC4 
>     28     22     101     25     4          5 
> #4392C6 #65AAD3 #69ACD5 #6CAED6 #77B4D8 #98C6DF 
>     14   6        8        4      3       9 
> plot(move, col=movers.den, pch=20,ylab="y coordinate movement (meters)",xlab="x coordinate movement (meters)")
> abline(h=0, v=0, col = "grey", lty=2)
> #legend??
> 
> Any help would be appreciated!
> Kate
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.