sm.density.compare():在单个图中显示多个密度估计
Posted
技术标签:
【中文标题】sm.density.compare():在单个图中显示多个密度估计【英文标题】:sm.density.compare(): displaying multiple density estimations in a single plot 【发布时间】:2016-10-24 20:58:42 【问题描述】:我正在尝试在 R 中叠加三个不同的密度图,以创建一个显示所有三条线的图(叠加)。我已经安装/加载了sm
包,但我尝试将它与我的数据一起使用,但无济于事。
我简单地使用density()
创建了三个单独的数据图并绘制了值。我的代码如下所示:
library(sm)
set.seed(0)
x <- rnorm(100, 0, 1)
y <- rnorm(126, 0.3, 1.2)
z <- rnorm(93, -0.5, 0.7)
dx <- density(x)
dy <- density(y)
dz <- density(z)
plot(dx)
plot(dy)
plot(dz)
但是当我尝试使用sm.density.compare()
覆盖图表时:
sm.density.compare(dx,dy,model="equal")
我收到一条错误消息:
sm.density.compare(dx, dy, model = "equal") 中的错误: sm.density.compare 只能处理一维数据跟踪:
有人知道我该如何解决这个问题吗?我研究了很多,但没有成功。我对 R 还很陌生,真的可以使用帮助。
【问题讨论】:
不!它的完美,它的工作!谢谢! 【参考方案1】:如果你想使用sm.density.compare()
,那么不要使用density()
。
sm.density.compare()
本身正在做密度估计。具体来说,它是对分组数据进行密度估计,以便您可以在同一张图上绘制不同组的密度。
这是你真正需要做的:
## three groups, each of length: length(x), length(y), length(z)
group.index <- rep(1:3, c(length(x), length(y), length(z)))
## collect data together and use sm.density.compare()
den <- sm.density.compare(c(x,y,z), group = group.index, model = "equal")
## plot will be generated automatically
使用model = "equal"
时,sm.density.compare()
有返回值。看看str(den)
:
List of 4
$ p : num 0
$ estimaate : num [1:3, 1:50] 2.37e-07 3.81e-06 6.06e-10 2.17e-06 2.26e-05 ...
$ eval.points: num [1:50] -4.12 -3.94 -3.76 -3.58 -3.4 ...
$ h : num 0.376
h
包含用于所有密度估计的带宽,eval.points
包含估计点,而estimaate
是密度估计值的矩阵。 (Adrian 这里有个错字,应该是“estimate”,而不是“estimaate”,哈哈)。
sm
包中的所有函数,以前缀 sm.
开头,接受可选参数 ...
,传递给 sm.options
。阅读?sm.options
,您会发现您对颜色显示、线型和线宽、带宽选择方法等拥有完全的控制权。
参考带只会添加到两组数据中。即,对于成对比较,sm.density.compare()
可以做得更多。例如:
den2 <- sm.density.compare(c(x,y), group = rep(1:2, c(length(x), length(y))),
model = "equal")
> str(den2)
List of 6
$ p : num 0.22
$ estimate : num [1:2, 1:50] 4.92e-06 2.70e-05 2.51e-05 1.00e-04 1.09e-04 ...
$ eval.points: num [1:50] -4.12 -3.94 -3.76 -3.58 -3.4 ...
$ upper : num [1:50] 0.00328 0.00373 0.00459 0.00614 0.00886 ...
$ lower : num [1:50] 0 0 0 0 0 ...
$ h : num 0.44
其中lower
和upper
给出了参考带/置信区域的界限。
如果您使用density()
,则不要使用sm.density.compare()
## set universal estimation range
xlim <- range(x, y, z)
dx <- density(x, from = xlim[1], to = xlim[2], n = 200)
dy <- density(y, from = xlim[1], to = xlim[2], n = 200)
dz <- density(z, from = xlim[1], to = xlim[2], n = 200)
在这种情况下,每个组的密度估计都是独立完成的。每个“密度”对象都是一个列表,例如:
> str(dx)
List of 7
$ x : num [1:200] -2.64 -2.61 -2.58 -2.55 -2.52 ...
$ y : num [1:200] 0.023 0.026 0.0291 0.0323 0.0356 ...
$ bw : num 0.31
$ n : int 100
$ call : language density.default(x = x, n = 200, from = xlim[1], to = xlim[2])
$ data.name: chr "x"
$ has.na : logi FALSE
- attr(*, "class")= chr "density"
x
是评估点,y
是估计密度,bw
是使用的带宽。您会看到,dx$bw
、dy$bw
和 dz$bw
是不同的,因为它们是独立估计的。但是,您可以在调用density()
时通过使用参数bw
手动指定通用bw
。见?density
,这里就不举例了。
现在,要覆盖这些密度图,您需要自己动手。
## set global plotting range
ylim <- range(dx$y, dy$y, dz$y)
## make plot
plot(dx$x, dx$y, col = 1, lwd = 2, type = "l", xlim = xlim, ylim = ylim)
lines(dy$x, dy$y, col = 2, lwd = 2)
lines(dz$x, dz$y, col = 3, lwd = 2)
【讨论】:
以上是关于sm.density.compare():在单个图中显示多个密度估计的主要内容,如果未能解决你的问题,请参考以下文章