阴影密度曲线下的区域，以标记最高密度间隔 (HDI)答案

【问题标题】：Shade an area under density curve, to mark the Highest Density Interval (HDI)阴影密度曲线下的区域，以标记最高密度间隔 (HDI)
【发布时间】：2020-04-15 17:24:25
【问题描述】：

我认为这应该很简单，但我迷路了，尽管网上有大量信息。

我的问题：我有一个数据点向量，我想为其绘制密度曲线，然后为曲线下的区域着色以表示最高密度间隔 (HDI)。当然，我正在尝试使用 ggplot2 包，特别是 qplot() 来实现这一点，因为我的数据是一个向量，而不是一个数据框。

可重现的示例

library(ggplot2)
library(HDInterval)

## create data vector
set.seed(789)
dat <- rnorm(1000)

## plot density curve with qplot and mark 95% hdi
qplot(dat, geom = "density")+ 
  geom_vline(aes(xintercept = c(hdi(dat))))

所以我明白了：

但我真正想要的是这样的：

有没有一种简单的方法可以通过ggplot2::qplot 实现这一目标？

【问题讨论】：

标签： r ggplot2 ggridges

【解决方案1】：

您可以使用 ggridges 软件包执行此操作。诀窍是我们可以将HDInterval::hdi 作为分位数函数提供给geom_density_ridges_gradient()，并且我们可以用它生成的“分位数”来填充。 “分位数”是下尾、中间和上尾的数字。

作为一般建议，我建议不要使用qplot()。这更有可能引起混乱，将向量放入 tibble 并不费力。

library(tidyverse)
library(HDInterval)
library(ggridges)
#> 
#> Attaching package: 'ggridges'
#> The following object is masked from 'package:ggplot2':
#> 
#>     scale_discrete_manual

## create data vector
set.seed(789)
dat <- rnorm(1000)

df <- tibble(dat)

## plot density curve with qplot and mark 95% hdi
ggplot(df, aes(x = dat, y = 0, fill = stat(quantile))) + 
  geom_density_ridges_gradient(quantile_lines = TRUE, quantile_fun = hdi, vline_linetype = 2) +
  scale_fill_manual(values = c("transparent", "lightblue", "transparent"), guide = "none")
#> Picking joint bandwidth of 0.227

^{由reprex package (v0.3.0) 于 2019 年 12 月 24 日创建}

scale_fill_manual() 中的颜色按三组的顺序排列，因此，例如，如果您只想遮蔽左尾，则应写为values = c("lightblue", "transparent", "transparent")。

【讨论】：

【解决方案2】：

当我阅读这篇文章时，我非常感谢您的回答，Wilke。但我想知道如何调整 hdi 的可信质量。最后我找到了解决方案！当我弄清楚分位数参数的来源（我通过点 [[2]] 访问它）时，它单击了触发器。我编写了以下函数（因为将值传递给 HDInterval::hdi 无法开箱即用）：

hdi_custWidth <- function(...) {
  dots <- list(...)
  quantiles <- dots[[2]]
  hdi_width <- quantiles[[length(quantiles)]] # uses the last entry if its a vector which should be the biggest one; better pass a single double < 1.0
  if (is.na(hdi_width)) hdi_width <- .89 # happens is quantiles = 1L
  message(paste0('HDI credible interval width = ', hdi_width))
  HDInterval::hdi(dots[[1]], credMass = hdi_width)
}

您可以使用它来更改上面帖子中的 repex：

library(tidyverse)
library(HDInterval)
library(ggridges)
#> 
#> Attaching package: 'ggridges'
#> The following object is masked from 'package:ggplot2':
#> 
#>     scale_discrete_manual

## create data vector
set.seed(789)
dat <- rnorm(1000)

df <- tibble(dat)

## plot density curve with qplot and mark 95% hdi
ggplot(df, aes(x = dat, y = 0, fill = stat(quantile))) + 
  geom_density_ridges_gradient(quantile_lines = TRUE, quantile_fun = hdi_custWidth, quantiles = .90, vline_linetype = 2) +
  scale_fill_manual(values = c("transparent", "lightblue", "transparent"), guide = "none")
#> Picking joint bandwidth of 0.227

当然，您可以在分位数参数中选择 0 到 1 之间的任何值（不仅仅是 .90）并获得相应的 hdi 可信质量。

【讨论】：