sm.density.compare(): отображение множественных оценок плотности в одном графике

Question

Я пытаюсь наложить три разных графика плотности в 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 может обрабатывать только трассировку 1-d данных:

Кто-нибудь знает, как я могу это исправить? Я исследовал довольно много, но без успеха. Я довольно новичок в R и действительно могу использовать помощь.

Answer 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 матрица, из оценочных значений плотности. (У Адриана есть опечатка, это должно быть «оценка», а не «оценка», LOL).

Все функции из sm пакета, начиная с префиксом sm. принимают необязательные аргументы ..., переданных sm.options. Прочтите на ?sm.options, и вы обнаружите, что у вас есть полный контроль над цветным дисплеем, типом линии и шириной линии, методом выбора полосы пропускания и т. Д.

Опорная полоса будет добавляться только к данным двух групп. I.e., для парного сравнения, 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 отличаются друг от друга, благодаря независимой оценке. Однако вы можете вручную указать универсальный bw при вызове density(), используя аргумент 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)