полоса пропускания плотности ядра в R

У меня есть два вектора: 1) ~ 1000 выборочных средних и 2) соответствующие ~ 1000 стандартных отклонений этих средних. Я хотел бы создать график плотности ядра этих данных, используя средние выборки в качестве наблюдений, по которым оценивается плотность, и стандартные отклонения каждого среднего в качестве ширины полосы для каждого наблюдения. Проблема в том, что плотность позволяет использовать только вектор длины 1 в качестве полосы пропускания. Например:

plot(density(means,bw=error))

возвращает следующие предупреждения:

1: In if (!is.finite(bw)) stop("non-finite 'bw'") :
  the condition has length > 1 and only the first element will be used
2: In if (bw <= 0) stop("'bw' is not positive.") :
  the condition has length > 1 and only the first element will be used
3: In if (!is.finite(from)) stop("non-finite 'from'") :
  the condition has length > 1 and only the first element will be used
4: In if (!is.finite(to)) stop("non-finite 'to'") :
  the condition has length > 1 and only the first element will be used

... и я получаю график, который использует ошибку первого элемента в списке как полосу пропускания для всех моих наблюдений.

Любые идеи о том, как я могу реализовать отдельную, определяемую пользователем полосу пропускания для каждого наблюдения, используемого для построения графика плотности ядра?


r bandwidth kernel-density
person user3525465    schedule 11.04.2014    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Не похоже, что density поддерживает такую ​​спецификацию полосы пропускания. Я полагаю, вы могли бы прокатиться самостоятельно

mydensity <- function(means, sds) {
  x <- seq(min(means - 3*sds), max(means + 3*sds), length.out=512)
  y <- sapply(x, function(v) mean(dnorm(v, means, sds)))
  cbind(x, y)
}

Это будет намного медленнее, чем реальная функция (которая, похоже, использует fft в вычислениях). Вот он в работе, с небольшой полосой пропускания слева и большой справа:

set.seed(144)
means <- runif(1000)
sds <- ifelse(means < 0.5, 0.001, 0.05)
plot(mydensity(means, sds))

введите описание изображения здесь

person josliber♦    schedule 12.04.2014