更直观的密度估计图形：ggdensity的改进
ggdensity是一个新的ggplot2扩展包，用于展示二维密度估计，使用的方法是基于最高密度区域（HDR）的密度估计方法。（什么是HDR？简单的说就是在指定概率所覆盖的样本空间所有可能的区域中，HDR具有可能的最小区域。）
1 2d密度估计的hdr密度图
1-1 geom_hdr( )
geom_hdr( )函数执行2D密度估计，计算并绘制得到的最高密度区域的填充图。geom_hdr( )很直观，展示的是给定概率水平下的最小区域。这就是所说的最高密度区域（HDR）。默认情况下，图中展示的概率为50%、80%、95%和99%。
这个函数参数众多，我们来重点关注其中的估计方法参数method。这个参数有四个取值："kde"，"histogram"，"freqpoly"或"mvnorm".在下面的例子里，使用随机生成的样本，我们对比四种估计方法的结果：
library(tidyverse)
library(ggdensity)
library(patchwork)
theme_set(theme_minimal())
df <- data.frame(x = rnorm(1000), y = rnorm(1000))
p<- ggplot(df, aes(x, y)) + coord_equal()
p1<-p + geom_hdr()
p2<-p + geom_hdr(method = "mvnorm")
p3<-p + geom_hdr(method = "freqpoly")
p4<-p + geom_hdr(method = "histogram")
(p1+p2)/(p3+p4)
图-1
geom_hdr( ) 对比geom_density_2d_filled( )
ggplot2中绘制两个连续变量联合分布的标准方法是geom_density_2d( )或者geom_density_2d_filled( )。geom_density_2d_filled( )绘制的等高线是估计二元密度的等距水平集合，也就是以等距的高度获得三维曲面的水平切片。所以，我们可以观察到密度高低的差异，与geom_hdr()相比进行解释要困难的多。
下面的例子使用一个模拟的例子对比两个函数的效果：
df <- data.frame(x = rnorm(1000), y = rnorm(1000))
p5 <- ggplot(df, aes(x, y)) + coord_equal()+
geom_density_2d_filled()
p6 <- ggplot(df, aes(x, y)) + coord_equal()+
geom_hdr()
p5+p6
图-2
表示区域概率的计算变量probs是geom_hdr( )使用底层stat函数创建的，可以使用after_stat( )按照ggplot2中对计算变量的标准方式来映射这个变量：
library(palmerpenguins)
ggplot(penguins, aes(flipper_length_mm, bill_length_mm, fill = species)) +
geom_hdr(
aes(alpha=0.8,fill = after_stat(probs)),
xlim = c(160, 240), ylim = c(30, 70)
)
图-3
1-2 geom_hdr_lines( )
geom_hdr_lines( )用最高密度的边界曲线代替对密度区域的填充：
ggplot(penguins, aes(flipper_length_mm, bill_length_mm, color = species)) +
geom_hdr_lines(xlim = c(160, 240), ylim = c(30, 70)) +
geom_point(size = 1)
图-4
1-3 geom_hdr_points( )
ggdensity包也提供了将密度图和散点图结合在一起的展示方式。这种展示最直接的方法就是将散点绘制在密度图上。另外一种方法是geom_hdr_points( )函数，直接在散点上绘制表示HDR概率的颜色，而不再画出区域：
p7<-ggplot(penguins, aes(flipper_length_mm, bill_length_mm, fill = species)) +
geom_hdr(xlim = c(160, 240), ylim = c(30, 70)) +
geom_point(shape = 21)
p8<-ggplot(penguins, aes(flipper_length_mm, bill_length_mm, fill = species)) +
geom_hdr_points(xlim = c(160, 240), ylim = c(30, 70))
p7+p8
图-5
对于数据观测很多，散点图存在明显重叠的情况，这种方法的优势是很明显的。diamonds数据集有5万多个观测，下面比较ggplot2的散点图geom_point( )和
geom_hdr_points( )的效果：
p_points <- ggplot(diamonds, aes(carat, price)) +
geom_point(alpha=0.2)
p_hdr_points <- ggplot(diamonds, aes(carat, price)) +
geom_hdr_points(alpha=0.2)
p_points + p_hdr_points
图-6
2 联合概率密度已知
有时候，我们知道联合概率密度函数，也可以绘制类似的HDR密度图，相应的函数有geom_hdr_fun( )、geom_hdr_lines_fun( )和geom_hdr_points_fun( )。这一类函数的特点是需要指定密度函数参数f。
下面两个例子，假设一个联合分布由两个独立的指数分布随机变量生成：
2-1geom_hdr_fun( )
f <- function(x, y) dexp(x) * dexp(y)
ggplot() +
geom_hdr_fun(fun = f, xlim = c(0, 10), ylim = c(0, 10))
图-7
2-2 geom_hdr_points_fun( )
f <- function(x, y) dexp(x) * dexp(y)
df <- data.frame(x = rexp(1000), y = rexp(1000))
ggplot(df, aes(x, y)) +
geom_hdr_fun(fun = f, xlim = c(0, 10), ylim = c(0, 10))+
geom_hdr_points_fun(fun = f, xlim = c(0, 10), ylim = c(0, 10))
图-8
3 地毯图geom_hdr_rug( )
地毯图是密度图的补充，可以展示联合分布的边际密度。geom_hdr_rug( )用于展示边际HDR，但是要注意，联合HDR并不是边际HDR的乘积：
ggplot(cars, aes(speed, dist)) +
geom_hdr() +
geom_point(color = "red") +
geom_hdr_rug()
图-9
