Как оптимизировать функцию, которая содержит интеграл в R или Python? - PullRequest
Новая
       88

Как оптимизировать функцию, которая содержит интеграл в R или Python?

0 голосов
/ 14 октября 2019

Впервые задаю вопрос здесь, но определенно выиграл от просмотра всех удивительных ответов на другие вопросы.

Итак, проблема в том, что я пытаюсь оптимизировать функцию через MLE в R (или Python тоже работает), а часть функции представляет собой интеграл, содержащий параметры, которые я оцениваю. Функция выглядит следующим образом:

Функция логарифмического правдоподобия

, где F (w) - CDF логарифмического нормального распределения с неизвестным средним значением и сигма кПо оценкам, лямбда и эта тоже. Все остальные параметры известны, включая w *. Я попытался оценить эту функцию, используя команду mle в R. Но она не смогла пройти из-за составной части. Может ли кто-нибудь помочь мне в этом? Мои коды ниже. Большое спасибо!

Определения f (w) и F (w)

Мои коды, которые не работали: 

  func <- function(mu, sigma_w){
F_w = (-.5*(((log(w_star)-mu)/sigma_w))^2)

} 

ll <- function(lambda,sigma_w,eta,mu) {
  R <- N*log(lambda) + N_u*log(1/sqrt(2*pi)*exp(integrate(func,lower = -Inf, upper = w_star, rel.tol=1e-5)$value)) - lambda*1/sqrt(2*pi)*exp(integrate(func,lower = -Inf, upper = w_star, rel.tol=1e-5)$value)*sum(T_ui) + N_u*log(eta) + sum(log(1/(sigma_w*w_i)*1/sqrt(2*pi)*exp(-.5*((log(w_i-mu)/sigma_w)^2)))) - N*log(eta+lambda*1/sqrt(2*pi)*exp(integrate(integrate(func,lower = -Inf, upper = w_star, rel.tol=1e-5)$value)))
  return(R)                           
}

fit <-mle(ll,start = list(lambda =0, sigma_w = 0, eta =0, mu = 0))

Сообщение об ошибке: 

Error in f(x, ...) : argument "sigma_w" is missing, with no default

1 Ответ

0 голосов
/ 14 октября 2019

Это прекрасный пример использования замыкания в R (функции, которые возвращают функции). Они служат именно для MLE, попробуйте этот код со своими значениями: 


library(stats4)

ll <- function(lambda, eta, sigma_w, mu, w_start, w_i, N, N_u )
{
  # fix the others params
  w_start <- w_start
  mu <- mu
  T_ui <- T_ui #your sample
  w_i <- w_i #another sample
  sigma_w <- sigma_w
  N <- N
  N_u <- N_u 
  #return another funtion
  function(lambda, eta)
  {
    F_star <- log( plnorm( q = w_star,
                           meanlog = mu, sdlog = sigma_w) ) 
    ll.scalar <- N * lambda + N_u*log( F_star) - lambda * F_star * sum(T_ui) +
      N_u* log(eta) + sum( log( dlnorm(x = w_i))) - N * log(eta + F_star)

    return(-ll.scalar) # negative for use in stats4::mle
  }
}

# your values example
w_star <- 1
mu <- 1
T_ui <- c(1,2,3) #your sample
w_i <- c(1,2,3)/10 #another sample
sigma_w <- 1
N <- 1 
N_u <- 1 
# You build the function with the closure
myLL <- ll(lambda , eta, sigma_w , mu, w_star , w_i , N, N_u  )


# now optimitation 

result <- mle(myLL, start= list(lambda=0,  eta=0))
result

```
...