О моделировании гамма-плотности

Я имитирую некоторые гамма-случайные числа

plot(density(rgamma(10000,8.1,rate=0.00510)),lwd=2,las=1,cex.axis=0.75,
main=expression(paste("Gamma Distribution with",' scale ',alpha," and  rate 
",beta)))

plot(density(rgamma(10000,2.1,rate=0.00110)),lwd=2,las=1,cex.axis=0.75,
 main=expression(paste("Gamma Distribution with",' scale ',alpha," and  rate 
",beta)))


plot(density(rgamma(10000,2.1,rate=110)),lwd=2,las=1,cex.axis=0.75,
 main=expression(paste("Gamma Distribution with",' scale ',alpha," and  rate 
",beta)))

Первый сюжет

Второй сюжет

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

Мне нужно смоделировать этот тип гаммы с некоторым хвостом, в среднем около 1200. Я выбирал случайные числа, чтобы получить эти значения, учитывая определение ожидания и дисперсии для гамма-распределения, но в первом случае я получаю отрицательные числа, Я не хочу этого. Во втором случае то же самое, но также на обоих графиках вероятность по оси y настолько мала, что я хотел бы увеличить эту вероятность, но я не знаю, как выбрать адекватные параметры для этого.

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


person user3483060    schedule 10.03.2018    source источник
comment
вы не получаете отрицательные числа.   -  person MichaelChirico    schedule 10.03.2018
comment
попробуйте использовать dgamma вместо rgamma.   -  person MichaelChirico    schedule 10.03.2018
comment
Попробуйте, например, x <- seq(0,5000, length = 10000); y <- dgamma(x, 8.1, .0051); plot(x,y, type = "l")   -  person hpesoj626    schedule 10.03.2018
comment
Высота кривой — это плотность, а не вероятность. Вероятность находится в области под кривой.   -  person pjs    schedule 10.03.2018


Ответы (1)


Вы не получаете отрицательных чисел в первом случае.

sum(rgamma(10000,8.1,rate=0.00510)<0)

# [1] 0

Что касается плотности, превышающей 1, то в этом нет ничего плохого. Это должна быть площадь под кривой, которая должна быть ровно 1.

Для построения графика распределения вы также можете использовать dgamma вместо rgamma, например:

x <- seq(0,.12, length = 10000)
y <- dgamma(x,2.1,rate=110)
plot(x,y, type = "l")

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

А вот сюжет с использованием ggplot2.

ggplot(data.frame(x,y), aes(x,y)) + 
  geom_line()

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

person hpesoj626    schedule 10.03.2018
comment
обратите внимание, что вы можете использовать from=0, чтобы предотвратить вычисление плотности для отрицательных значений, но вычисление плотности ядра вблизи жестких границ сложно... - person Ben Bolker; 11.03.2018
comment
Спасибо, @BenBolker. Я помню, как комментировал сразу после того, как отредактировал свой ответ, но почему-то комментарий отсутствует. - person hpesoj626; 12.03.2018