如何将等高线添加到 R 中的密度散点图？

Question

我有一个代码可以生成这样的密度散点图（在这里找到How to add an ellipse to a scatter plot colored by density in R?）：

## Data in a data.frame
x1 <- rnorm(n=1000, sd=2)
x2 <- x1*1.2 + rnorm(n=1000, sd=2)
df <- data.frame(x1,x2)

## Use densCols() output to get density at each point
x <- densCols(x1,x2, colramp=colorRampPalette(c("black", "white")))
df$dens <- col2rgb(x)[1,] + 1L

## Map densities to colors
cols <-  colorRampPalette(c("#FF3100", "#FF9400", "#FCFF00", 
                            "#45FE4F", "#00FEFF", "#000099"))(6)
df$col <- ifelse(df$dens >= 250, cols[1], ifelse(df$dens >= 200, cols[2], ifelse(df$dens >= 150, cols[3], ifelse(df$dens >= 100, cols[4], ifelse(df$dens >= 50, cols[5], cols[6])))))

## Plot it, reordering rows so that densest points are plotted on top
plot(x2~x1, data=df[order(df$dens),], pch=20, col=col, cex=1)

我想添加与下图完全相同的轮廓（在此处找到 https://fr.mathworks.com/matlabcentral/fileexchange/8577-scatplot）：

我必须在没有 ggplot2 的情况下这样做。

请问有人有意见吗？

Answer 1

这些等高线显示给定级别的数据估计密度。您可以使用MASS::kde2d() 计算密度估计值，然后使用contour() 绘制：

biv_kde <- MASS::kde2d(x1, x2)
contour(biv_kde, add = T)

编辑：

这是使用KernSmooth 获取彩色密度轮廓的另一种方法：

（请注意，我在生成数据之前使用了set.seed(123)。）

# estimate bandwidths
h <- c(dpik(df$x1), dpik(df$x2))

# obtain density estimatte
f1 <- bkde2D(df[, 1:2], bandwidth = h, gridsize = c(200, 200))

# setup vector of density levels to obtain contours for
contour_levels <- pretty(f1$fhat, 9)

# setup color palette
crp <- colorRampPalette(rev(c("#FF3100", "#FF9400", "#FCFF00", 
                              "#45FE4F", "#00FEFF", "#000099")))

# density based colors 
df$col <- densCols(x1, x2, bandwidth = h, colramp = crp)

# plot
plot(x2 ~ x1, data = df, pch = 20, col = col, cex = 0.5)
contour(f1$x1, f1$x2, f1$fhat, levels = contour_levels, col = crp(9),
        lwd = 2, add = T)