diffdensity.plot: Plot the difference of densities between two distributions 1...
In xrobin/xavamess: Useful Miscellaneous Functions

Description Usage Arguments Value See Also Examples

Plot the difference of densities between two distributions 1 and 2

diffdensity.plot(
  x1,
  x2,
  y1,
  y2,
  bw = bw.nrd,
  n = 100,
  ztransform = sqrt_transform,
  zlim = NULL,
  zlim.mult = 1,
  relative = FALSE,
  plot.d1 = FALSE,
  plot.d2 = FALSE,
  plot.diff = TRUE,
  draw.contour = TRUE,
  contour.levels = NULL,
  xlim = extendrange(range(c(x1, x2)), f = extendrange.x),
  ylim = extendrange(range(c(y1, y2)), f = extendrange.y),
  extendrange = 0,
  extendrange.x = extendrange,
  extendrange.y = extendrange,
  col1 = "blue",
  col2 = "red",
  bg.col = "white",
  ncol = 101,
  add = FALSE,
  ...
)

`x1, x2, y1, y2`	the x and y values for distributions 1 and 2.
`bw`	a bandwidth selection function such as `bw.nrd`; or a numeric vector for `kde2d`.
`n`	the number of grid points.
`ztransform`	transformation of the z axis of the 2D density. The default `sqrt_transform` will take square root the differences in `z` to increase contrast. Use `I` in order to not transform the axis.
`zlim`	overrides the limits in the z direction, setting the value where either `col1` or `col2` is displayed in full brightness.
`zlim.mult`	a multiplication factor for zlim, especially useful in `relative = TRUE` mode.
`relative`	whether to adjust `zlim` to match the highest density in distribution 1 and 2 (`TRUE`) or to the highest difference (`FALSE`). Ignored if `zlim` is provided.
`plot.d1, plot.d2`	whether to plot the 2d densities of distributions 1 and 2, respectively. Disabled by default.
`plot.diff`	whether to plot the differential densities. Enabled by default.
`draw.contour`	whether to overlay a `contour` on the diff plot. Ignored if `plot.diff = FALSE`.
`contour.levels`	levels at which to draw contour lines, see the `levels` argument to `contour`
`xlim, ylim`	limits of the plot on x and y.
`extendrange, extendrange.x, extendrange.y`	extend `xlim` and `ylim` by a this fraction. Ignored if `xlim` and `ylim` are provided.
`col1, col2`	colors to represent the two distributions.
`bg.col`	the background color.
`ncol`	the number of colors of the palette. Should be an odd number to avoid weird effects around 0
`add`	if TRUE, add to the current plot instead of creating a new one.
`...`	further arguments from and two other methods, in particular graphical parameters for `image` such as `main`, `xlab` or `ylab`.

Invisibly, a list like kde2d with z the difference in densities, and the following additional “z” elements:

d1.z, d2.zDensities of distributions 1 and 2
sum.zThe sum of the density distributions 1 and 2

The Stack Overflow answer that inspired this function

# Get some gaussian distributions
x1 <- rnorm(10000)
y1 <- rnorm(10000)
x2 <- rnorm(10000, mean = 1, sd = 0.5)
y2 <- rnorm(10000, mean = 0.5, sd = 2)
# First shot at the plot
diffdensity.plot(x1, x2, y1, y2)

# Make it smoother
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10)

# Improve contour
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
		contour.levels = c(0.01, 0.1, 0.2, 0.3, 0.4))

# Change the colors
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
		col1 = "green", col2 = "yellow", bg.col = "black")

# Also plot d1 and d2
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
plot.d1 = TRUE, plot.d2 = TRUE, relative = TRUE)

# Increase resolution
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
		n = 500)

# Try an assymetric distribution
y1 <- rexp(10000, rate = 1)
y2 <- rexp(10000, rate = 1.5)
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10)

# A bit off the density is clipped in y, so extend the limits
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
                 extendrange.y = .1)
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
                 extendrange = .1, extendrange.x = 0)
# Extend in x as well
diffdensity.plot(x1, x2, y1, y2, bw = function(x) bw.nrd(x) * 10,
                 extendrange = .1)