Оценка максимального правдоподобия для биномиальных данных - PullRequest
Новая
       6

Оценка максимального правдоподобия для биномиальных данных

0 голосов
/ 05 декабря 2018

Мне нужно найти оценку максимального правдоподобия для вектора биномиальных данных.

один такой: 

binvec <- rbinom(1000, 1, 0.5)

Я попытался сначала создать функцию, а затем оптимизировать ее с помощью optim(). 

llbinom <- function(theta, x) {return(sum(dbinom(x=x,
               size = theta[1],
               prob = theta[2],log=TRUE)))}
mybin <- optim(c(0.5,0.5),fn=llbinom,x=binvec)
mybin

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

Вот исходный код, который я получил: 

ll <- function(theta,x) {return(-sum(dnorm(x=x,
   mean=theta[1],sd=theta[2],log=TRUE)))}
mle <- optim(c(5,3),fn=ll,x=binvec)

1 Ответ

0 голосов
/ 05 декабря 2018

Несколько проблем здесь.

  • похоже, что вы пропускаете отрицательный знак (optim() по умолчанию минимизирует , если вы не установите параметр управления fnscale=-1, так что вынеобходимо определить отрицательная функция логарифмического правдоподобия)
  • параметр size должен быть целым числом
  • это необычно и технически сложно, чтобы оценить size параметр из данных (это часто делается с использованием моделей N-смеси , если вы хотите ознакомиться с методикой: например, посмотрите пакет unmarked);обычно число испытаний считается известным.Поэтому я бы попробовал 
llbinom <- function(theta, x) {return(
   -sum(dbinom(x=x,
               size = 1,
               prob = theta[1],log=TRUE)))}
mybin <- optim(c(0.5),fn=llbinom,x=binvec)

Есть много причин сделать это трудным путем (численно);если вам действительно нужно только найти MLE вероятности одной биномиальной выборки x (независимые наблюдения с той же вероятностью успеха из s испытаний), аналитическое решение будет sum(x)/sum(s) ...

...