Функция R Optim (): усеченная логарифмическая нормальная оценка максимального правдоподобия (решение для mu и sd) - PullRequest
0 голосов
/ 17 марта 2020

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

log.lklh.lnorm <- function(x, mu, sd, input_min, input_max){-sum(log(dlnorm(x, meanlog = mu, sdlog = sd, log = FALSE)/((plnorm(input_max, meanlog = mu, sdlog = sd, lower.tail = TRUE, log.p = FALSE)) - (plnorm(input_min, meanlog = mu, sdlog = sd, lower.tail = TRUE, log.p = FALSE)))))}


optim(par = c(0,1,1,100000), log.lklh.lnorm)

Я использовал Excel Solver, чтобы решить (т.е. найти максимальное значение этой суммы) для mu и sd (с учетом ввода max и min). Тем не менее, я не могу повторить это в R. Я пробовал различные версии приведенного выше кода, в том числе:

  • Установите для input min и input max определенные значения, такие как 1 и 100 000 соответственно.
  • Порядок вызовов Optim (), то есть порядок ввода x, par и fn)

Любая помощь очень ценится, спасибо!

Дейв

1 Ответ

1 голос
/ 19 марта 2020

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

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

log.lklh.lnorm <- function(x, mu, sd, input_min, input_max){-sum(log(dlnorm(x, meanlog = mu, sdlog = sd, log = FALSE)/((plnorm(input_max, meanlog = mu, sdlog = sd, lower.tail = TRUE, log.p = FALSE)) - (plnorm(input_min, meanlog = mu, sdlog = sd, lower.tail = TRUE, log.p = FALSE)))))}


f <- function(y) {
  log.lklh.lnorm(x,y[1],y[2],1,100000)    
}

optim(par = c(3,1), f)

Примечания в ответ на дополнительные вопросы в комментариях:

  1. Какова функция par = c(3,1)? Это имеет две функции. (1) Это начальная точка (предположение). (2) Определяет количество параметров для оптимизации. В данном случае это 2. Следует отметить, что вы можете (и должны) рассчитывать очень хорошие начальные значения. Вы можете получить очень хорошие оценки для μ и σ, просто рассчитав среднее значение и стандартное отклонение выборки. (В статистике это иногда называется method of moments). Это значительно облегчит задачу optim. Поэтому вместо par=c(3,1) мы действительно должны сделать: par = c(mean(x),sd(x)).

  2. Что делает y[1]/y[2]? Функция optim организует все параметры для оптимизации, как одиночный вектор . Таким образом, я представлял (μ, σ) одним вектором y (длиной 2). В функции f этот массив распаковывается и передается log.lklh.lnorm в качестве отдельных аргументов.

...