Найти интервалы плотности вероятности

9

У меня есть вектор

x <- c(1,2,3,4,5,5,5,6,6,6,6,
       7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
       7,7,7,7,7,7,7,7,8,8,8,8,9,9,9,10)

(мой фактический вектор имеет длину> 10000), и я хотел бы найти интервалы, в которых лежит 90% плотности. Является ли quantile(x, probs=c(0.05,0.95), type=5)наиболее подходящим или есть другой способ?

ECII
источник
Ваш вопрос немного неясен насчет "интервалов, где ..." - может быть несколько интервалов. Вас интересуют только внутренние 90%, то есть симметричная обрезка с каждой стороны? В конце концов, от минимального до 90% иль захватывает 90% данных, аналогично для 10% ил до максимального значения.
Итератор
Вы ищете кратчайший интервал, симметричный интервал (равная вероятность на каждом конце) или что-то еще?
Glen_b

Ответы:

19

Как указано выше, существует много разных способов определения интервала, который включает 90% плотности. Еще не указанным является самый высокий [задний] интервал плотности ( википедия ), который определяется как «кратчайший интервал, для которого разница в значениях эмпирической функции накопленной плотности конечных точек является номинальной вероятностью».

library(coda)
HPDinterval(as.mcmc(x), prob=0.9)
Бен Болкер
источник
3

Это, безусловно, кажется самым простым подходом. Функция довольно быстрая. Я использую его все время на выборках, которые в сотни раз больше той, которую вы используете, и стабильность оценок должна быть хорошей для вашего размера выборки.

В других пакетах есть функции, которые предоставляют более полные наборы описательной статистики. Я использую один Hmisc::describe, но есть несколько других пакетов с describeфункциями.

Dwin
источник
3

Ваш путь кажется разумным, особенно с дискретными данными в примере,

quantile(x,probs=c(0.05,0.95), type=5)
 5% 95% 
2.8 9.0

но другой способ - использовать вычисленное ядро ​​плотности:

dx <- density(x)
dn <- cumsum(dx$y)/sum(dx$y)
li <- which(dn>=0.05)[1]
ui <- which(dn>=0.95)[1]
dx$x[c(li,ui)]
[1] 2.787912 9.163246
Джеймс
источник
-1

Да. :-). Вы можете найти вывод stats::densityболее полезным.

Карл Виттофт
источник