如何在 R 中叠加密度图？

Posted 2023-02-16

【中文标题】如何在 R 中叠加密度图？

【英文标题】：How to overlay density plots in R?

【发布时间】：2011-10-19 19:53:21

【问题描述】：

我想用 R 在同一设备上叠加 2 个密度图。我该怎么做？我搜索了网络，但没有找到任何明显的解决方案。

我的想法是从文本文件（列）中读取数据，然后使用

plot(density(MyData$Column1))
plot(density(MyData$Column2), add=T)

或者本着这种精神的东西。

【问题讨论】：

对于ggplot2 家族，现在有一个包“ggridges”可以做到这一点。

【参考方案1】：

每当出现轴限制不匹配的问题时，base 图形中的正确工具是使用matplot。关键是利用from 和to 参数到density.default。这有点骇人听闻，但很容易自己动手：

set.seed(102349)
x1 = rnorm(1000, mean = 5, sd = 3)
x2 = rnorm(5000, mean = 2, sd = 8)

xrng = range(x1, x2)

#force the x values at which density is
#  evaluated to be the same between 'density'
#  calls by specifying 'from' and 'to'
#  (and possibly 'n', if you'd like)
kde1 = density(x1, from = xrng[1L], to = xrng[2L])
kde2 = density(x2, from = xrng[1L], to = xrng[2L])

matplot(kde1$x, cbind(kde1$y, kde2$y))

根据需要添加花里胡哨（matplot 接受所有标准的plot/par 参数，例如lty、type、col、lwd，...）。

【参考方案2】：

您可以使用ggjoy 包。假设我们有三个不同的beta 分布，例如：

set.seed(5)
b1<-data.frame(Variant= "Variant 1", Values = rbeta(1000, 101, 1001))
b2<-data.frame(Variant= "Variant 2", Values = rbeta(1000, 111, 1011))
b3<-data.frame(Variant= "Variant 3", Values = rbeta(1000, 11, 101))


df<-rbind(b1,b2,b3)

你可以得到如下三种不同的分布：

library(tidyverse)
library(ggjoy)


ggplot(df, aes(x=Values, y=Variant))+
    geom_joy(scale = 2, alpha=0.5) +
    scale_y_discrete(expand=c(0.01, 0)) +
    scale_x_continuous(expand=c(0.01, 0)) +
    theme_joy()

【参考方案3】：

这就是我在 base 中的做法（实际上在第一个答案 cmets 中提到了，但我会在这里展示完整的代码，包括图例，因为我还不能评论......）

首先，您需要从密度图中获取有关 y 轴最大值的信息。所以你需要先分别计算密度

dta_A <- density(VarA, na.rm = TRUE)
dta_B <- density(VarB, na.rm = TRUE)

然后根据第一个答案绘制它们，并为刚刚得到的 y 轴定义最小值和最大值。 （我将最小值设置为 0）

plot(dta_A, col = "blue", main = "2 densities on one plot"), 
     ylim = c(0, max(dta_A$y,dta_B$y)))  
lines(dta_B, col = "red")

然后在右上角添加一个图例

legend("topright", c("VarA","VarB"), lty = c(1,1), col = c("blue","red"))

【参考方案4】：

第二个使用lines：

plot(density(MyData$Column1))
lines(density(MyData$Column2))

不过，请确保第一个图的范围是合适的。

【讨论】：

+1 当两个密度具有不同的范围并且第二条曲线不适合绘图限制时，您可能需要稍微复杂一些的东西。然后您可以在绘图之前计算密度，并使用range(dens1$y, dens2$y) 计算适当的ylim，其中dens1 和dens2 是包含两个密度估计对象的对象。在对plot() 的调用中使用此ylim。

您可能还想区分这两行。设置线条宽度 (lwd)、线条类型 (lty) 或线条颜色 (col) 应该会有所帮助。此时，您还可以考虑添加一个图例，使用 legend()

@Gavin 如果 OP 正在从文件中读取，我将构建一个复杂的函数来读取数据（sapply，lapply），查找所有数据集的范围，将默认范围设置为最大值所有的范围，然后绘制（线）密度。

【参考方案5】：

添加处理 y 轴限制的基本图形版本，添加颜色并适用于任意数量的列：

如果我们有一个数据集：

myData <- data.frame(std.nromal=rnorm(1000, m=0, sd=1),
                     wide.normal=rnorm(1000, m=0, sd=2),
                     exponent=rexp(1000, rate=1),
                     uniform=runif(1000, min=-3, max=3)
                     )

然后绘制密度：

dens <- apply(myData, 2, density)

plot(NA, xlim=range(sapply(dens, "[", "x")), ylim=range(sapply(dens, "[", "y")))
mapply(lines, dens, col=1:length(dens))

legend("topright", legend=names(dens), fill=1:length(dens))

这给出了：

【讨论】：

我喜欢这个例子，但如果你有包含 NA 值的数据列，它就不起作用。我不确定如何修改代码，但这会很有用

@daisy 将此行 dens <- apply(myData, 2, density) 更改为 dens <- apply(myData, 2, density, na.rm=TRUE)，它应该可以工作。

【参考方案6】：

我拿上面的格子例子做了一个漂亮的函数。可能有更好的方法通过熔体/铸造进行重塑。 （如果您看到改进，请发表评论或编辑。）

multi.density.plot=function(data,main=paste(names(data),collapse = ' vs '),...)
  ##combines multiple density plots together when given a list
  df=data.frame();
  for(n in names(data))
    idf=data.frame(x=data[[n]],label=rep(n,length(data[[n]])))
    df=rbind(df,idf)
  
  densityplot(~x,data=df,groups = label,plot.points = F, ref = T, auto.key = list(space = "right"),main=main,...)

示例用法：

multi.density.plot(list(BN1=bn1$V1,BN2=bn2$V1),main='BN1 vs BN2')

multi.density.plot(list(BN1=bn1$V1,BN2=bn2$V1))

【参考方案7】：

只是为了提供一个完整的集合，这是使用lattice的Chase答案的一个版本：

dat <- data.frame(dens = c(rnorm(100), rnorm(100, 10, 5))
                   , lines = rep(c("a", "b"), each = 100))

densityplot(~dens,data=dat,groups = lines,
            plot.points = FALSE, ref = TRUE, 
            auto.key = list(space = "right"))

产生这样的情节：

【讨论】：

不创建新的data.frame:densityplot(~rnorm(100)+rnorm(100, 10, 5), plot.points=FALSE, ref=TRUE, auto.key = list(space = "right"))。或者对于 OP 数据densityplot(~Column1+Column2, data=myData).

【参考方案8】：

ggplot2 是另一个图形包，它以一种非常巧妙的方式处理 Gavin 提到的范围问题。它还可以处理自动生成适当的图例，并且在我看来，开箱即用的感觉通常更精致，手动操作更少。

library(ggplot2)

#Sample data
dat <- data.frame(dens = c(rnorm(100), rnorm(100, 10, 5))
                   , lines = rep(c("a", "b"), each = 100))
#Plot.
ggplot(dat, aes(x = dens, fill = lines)) + geom_density(alpha = 0.5)

【讨论】：

OP 的 data.frame 需要先熔成长格式：ggplot (melt (MyData), mapping = aes (fill = variable, x = value)) + geom_density (alpha = .5)

不错的情节。什么是“dat2”...？什么是“融化”（未找到命令）？

@ErikAronesty - 你猜在这一点上和我的一样好，我在两年前回答了这个问题！我推测我的环境中有另一个名为dat 的对象，因此将其命名为dat2...我提供的模拟数据与宣传的一样。 melt() 命令来自包 reshape2。早在 2011 年，reshape2 会在 ggplot2 加载时自动加载，但现在情况已不再如此，因此您需要单独执行 library(reshape2)。