在 ggplot2 中使用分面绘制分布的分位数。

Question

我目前正在使用ggplot绘制多个回归模型的一阶差分分布。为了便于解释差异，我想标记每个分布的2.5％和97.5％百分位数。由于我将要做很多图，并且数据分组在两个维度（模型和类型）中，我想在ggplot环境中定义和绘制相应的百分位数。使用facet将分布绘制到我想要的位置，但不包括百分位数。当然，我可以更加手动地进行操作，但我希望能找到一个解决方案，在这个解决方案中我仍然能够使用facet_grid,因为这样可以减少与拟合不同图形一起的麻烦。

下面是使用模拟数据的示例：

df.example <- data.frame(model = rep(c("a", "b"), length.out = 500), 
                      type = rep(c("t1", "t2", "t2", "t1"), 
                      length.outh = 250), value = rnorm(1000))

 ggplot(df.example, aes(x = value)) +
 facet_grid(type ~ model) +
 geom_density(aes(fill = model, colour = model))

我尝试过两种方法来添加分位数。第一种方法会产生错误信息：

 ggplot(df.example, aes(x = value)) +
 facet_grid(. ~ model) +
 geom_density(aes(fill = model, colour = model)) +
 geom_vline(aes(x = value), xintercept = quantile(value, probs = c(.025, .975)))

Error in quantile(value, probs = c(0.025, 0.975)) : object 'value' not found

第二个选项可以为整个变量获取分位数，而不是子密度的分位数。也就是说，所有四个密度的绘制分位数是相同的。

 ggplot(df.example, aes(x = value)) +
 facet_grid(type ~ model) +
 geom_density(aes(fill = model, colour = model)) +
 geom_vline(xintercept = quantile(df.example$value, probs = c(.025, .975)))

我想知道在ggplot2环境中是否有一种方法来绘制每个子组的特定分位数？ 感谢任何意见。

- chrstnsn

Answer 1

使用plyr（或dplyr、data.table）预计算这些值...

set.seed(1)
# ...

df.q <- ddply(df.example, .(model, type),
              summarize, q=quantile(value, c(.025, .975)))    
p + geom_vline(aes(xintercept=q), data=df.q)

- rcs

---

非常好的建议。谢谢！ - chrstnsn

Answer 2

您可以预先计算分位数。

使用您的示例数据：

library (dplyr)
d2 <- df.example %>%
  group_by(model, type) %>%
  summarize(lower = quantile(value, probs = .025),
            upper = quantile(value, probs = .975))

然后像这样绘制：

ggplot(df.example, aes(x = value)) +
  facet_grid(type ~ model) +
  geom_density(aes(fill = model, colour = model)) +
  geom_vline(data = d2, aes(xintercept = lower)) +
  geom_vline(data = d2, aes(xintercept = upper))

- Axeman

---

非常好的建议。谢谢！ - chrstnsn

Answer 3

现在，可以使用orientation选项与stat_summary()一起使用，无需预先计算即可实现相同的结果。

为每个面板定义一个虚拟的y值，以及使用orientation = "y"来对观测值进行分组。然后，在stat_summary()中使用自定义fun来计算每个面板所需的分位数向量。要将结果显示为垂直线，请在美学规范中指定geom = "vline"和其所需的xintercept，并使用after_stat(x)从计算出的x值中获取结果，现在使用fun计算得出的结果。

library(ggplot2)

set.seed(1)

df.example <- data.frame(
  model = rep(c("a", "b"), length.out = 500),
  type = rep(c("t1", "t2", "t2", "t1"),
    length.outh = 250
  ), value = rnorm(1000)
)

ggplot(df.example, aes(x = value)) +
  facet_grid(type ~ model) +
  geom_density(aes(fill = model, colour = model)) +
  stat_summary(
    geom = "vline",
    orientation = "y",
    # y is a required aesthetic, so use a dummy value
    aes(y = 1, xintercept = after_stat(x)),
    fun = function(x) {
      quantile(x, probs = c(0.025, 0.975))
    }
  )

- Mikko Marttila

---

如果您使用fun而不是fun.data，则可以使代码稍微简单一些：stat_summary(aes(y = 1, xintercept = after_stat(x)), fun = function(x) quantile(x, probs = c(0.025, 0.975)), geom = "vline", orientation = "y") - Axeman

或者 stat_summary(aes(y = 1, xintercept = after_stat(x)), fun = quantile, fun.args = list(probs = c(0.025, 0.975)), geom = "vline", orientation = "y")。 - Axeman

@Axeman 啊是啊，谢谢你指出来。那样会更简单。现在另一个想法是保留 fun.data 但是改变数据框架使其拥有一个名为 yintercept 的列，并且完全省掉 aes() 中的 xintercept = stat(x)。 - Mikko Marttila

话虽如此，我认为“fun”方法肯定更容易记住。我已经进行了编辑，仍然偏爱匿名函数，因为我个人觉得它更易读。再次感谢！ - Mikko Marttila

Answer 4

好问题。同样的问题的更一般版本是：在使用 facets 时如何在子集数据集上调用函数？这似乎是一个非常有用的功能，所以我搜索了一下，但没有找到任何相关信息。

已经给出的答案非常好。另一个选项是使用 multiplot() 作为手动分面的一种方式。

- CoderGuy123

---

我同意。这两个解决方案都非常简洁，但正如你所指出的那样，并没有真正解决我在问题中提出的问题。这确实将是ggplot一个非常有趣的功能。 - chrstnsn