Подбор логнормального распределения или распределения Пуассона

У меня есть вектор из 1096 чисел, среднесуточная концентрация NOx, измеренная за 3 года на измерительной станции. Вы можете наблюдать тип распределения на изображении:

Гистограмма концентрации NOx

Я использовал эти команды для построения гистограммы:

NOxV<-scan("NOx_Vt15-17.txt")
hist.NOxVt<-hist(NOxV, plot = FALSE, breaks = 24) 
plot(hist.NOxVt, xlab = "[NOx]", ylab = "Frequenze assolute", main = "Istogramma freq. ass. NOx 15-17 Viterbo")
points(hist.NOxVt$mids, hist.NOxVt$counts, col= "red")

Мой профессор предложил мне подогнать гистограмму к распределению Пуассона — обращая внимание на переход: дискретное -> непрерывное (я не знаю, что это значит) — или к «логарифмически нормальному» распределению.

Я попытался выполнить подгонку Пуассона с некоторыми командными строками, которые она дала нам на уроке, но R выдал мне ошибку после выполнения последней строки кода из следующего:

  my_poisson = function(params, x){
      exp(-params)*params^x/factorial(x)
  }

  y<-hist.NOxVt$counts/1096;
  x<-hist.NOxVt$mids;
  z <- nls( y ~ exp(-a)*a^x/factorial(x), start=list(a=1) )

Ошибка в numericDeriv(form[[3L]], имена(ind), env): пропущенное значение или бесконечность, полученная при оценке модели. Кроме того: было 50 или более предупреждений (используйте warnings(), чтобы увидеть первые 50)"

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

Я был бы признателен за любые предложения или примеры того, как сделать логнормальную подгонку и/или подгонку Пуассона.


r model-fitting poisson data-fitting normal-distribution
person Silvia    schedule 11.03.2018    source источник

Ответы (1)


arrow_upward
5
arrow_downward

В пакете MASS, который поставляется с R, есть встроенная функция fitdistr:

Создание примера данных для просмотра (параметры зрительного восприятия, чтобы получить что-то похожее на вашу картинку):

set.seed(101)
z <- rlnorm(1096,meanlog=4.5,sdlog=0.8)

Подгонка (по статистическим соображениям я бы не рекомендовал подгонку Пуассона - может быть возможно адаптировать дискретное распределение, такое как Пуассона (или, лучше, отрицательное биномиальное), чтобы соответствовать таким непрерывным данным, но логарифмически-нормальное или гамма-распределения являются более естественным выбором. .

library(MASS)
f1 <- fitdistr(z,"lognormal")
f2 <- fitdistr(z,"Gamma")

Объекты f1 и f2 при печати дают расчетные коэффициенты (meanlog и sdlog для логарифмически-нормального, shape и rate для гаммы) и стандартные ошибки коэффициентов.

Нарисуйте картинку (по шкале плотности, а не по шкале счета): красный — логарифмическая норма, синий — гамма (в данном случае логарифмическая норма подходит лучше, потому что именно так я сгенерировал «данные» в первую очередь). [Материал with(as.list(coef(...)) — это некоторая причудливость R, позволяющая использовать имена коэффициентов (meanlog, sdlog и т. д.) в последующем коде R.]

hist(z,col="gray",breaks=50,freq=FALSE)
with(as.list(coef(f1)),
     curve(dlnorm(x,meanlog,sdlog),
           add=TRUE,col="red",lwd=2))
with(as.list(coef(f2)),
     curve(dgamma(x,shape=shape,rate=rate),
           add=TRUE,col="blue",lwd=2))

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

person Ben Bolker    schedule 11.03.2018