优化最大似然函数

Question

我从正态分布生成了一个包含 20 个随机点的数据集，创建了对应于这 20 个点的最大似然函数，然后尝试优化这个函数以找出均值 (mu) 和标准差 (sigma) ).

首先，我生成了随机数据：

y <- rnorm(20,5,5)

然后，我定义了最大似然函数：

my_function <- function(x) {
 
  n = 20
  
  a = -n/2*log(2*pi)
  b = -n/2*log(x[2]^2)
  c = -1/(2*x[2]^2)
  d = sum ((y-x[1])^2)     
  return( a + b + c* d)
  
}

然后，我尝试用三种不同的方式优化这个函数 - 但所有这些方式都失败了：

方法一：

#load all libraries
library(pracma)
library(stats)
library(GA)

    steep_descent(c(1, 1), my_function,  maxiter = 10000)
   
$xmin
[1] NA

$fmin
[1] NA

$niter
[1] 10001

Warning message:
In steep_descent(c(1, 1), my_function, maxiter = 10000) :
  Maximum number of iterations reached -- not converged.

方法二：

 optim(c(1,1), my_function)

$par
[1] -8.487500e+00  1.776357e-16

$value
[1] -5.559631e+34

$counts
function gradient 
     501       NA 

$convergence
[1] 1

$message 
NULL

方法三：

# I redefined the function for the GA Library
my_function <- function(x_1, x_2) {
 
  n = 20
  
  a = -n/2*log(2*pi)
  b = -n/2*log(x_2^2)
  c = -1/(2*x_2^2)
  d = sum ((y-x_1)^2)     
  return( a + b + c* d)
  
}


GA <- ga(type = "real-valued", 
         fitness = function(x)  my_function(x[1], x[2]),
         lower = c(-20,20), upper = c(-20,20), 
         popSize = 50, maxiter = 1000, run = 100)

GA@solution
      x1 x2
[1,] -20 20

在所有 3 种方法中，从未找到正确的正确答案 - 事实上，3 种方法中的 none 找到了接近正确答案的任何东西。

根据统计理论，正确答案（基于此问题随机生成的数据）应该是：

> mu = mean(y)
> sigma = sd(y)

> mu
[1] 3.937513

> sigma
[1] 4.707227

我做错了什么，我该如何解决这个问题（即优化方法 return 答案更接近正确答案）？

Answer 1

optim和steep_descent都执行最小化（你想要最大化）。在 GA 中，您的参数是固定的（例如 lower = upper = 20）。调整您的代码：

my_function <- function(param, vec) {
  
  - length(vec)/2 * log(param[[2]]) - 1/(2 * param[[2]]) * sum((vec - param[[1]])^2)
}

优化

res1 <- optim(
  method = "L-BFGS-B",
  par = c(0, 1),
  fn = my_function,
  vec = y,
  lower = c(-20, 0.01),
  upper = c(20, 40),
  control = list(fnscale = -1)
)

fnscale 的负值将执行最大化。

GA

res2 <- GA::ga(
  type = "real-valued",
  popSize = 100,
  fitnes = my_function,
  lower = c(-20, 0),
  upper = c(20, 40),
  vec = y,
  maxiter = 1000,
  run = 100,
  optim = TRUE
)

为了使用本地搜索方法，最好放optim=TRUE（包只支持real-valued问题）

DEoptim

你也可以使用差分进化（它也做最小化）：

res3 <- DEoptim::DEoptim(
  fn = function(param, vec) -my_function(param, vec),
  lower = c(-20, 0),
  upper = c(20, 40),
  vec = y,
  control=list(NP = 100, itermax = 1000))
)

结果：

> res1$par
[1]  1.935102 18.683264
> res2@solution[1,]
       x1        x2 
 1.935101 18.683222 
> res3$optim$bestmem
     par1      par2 
 1.935101 18.683252 
>

（为了比较结果需要种子）

优化最大似然函数

Optimizing Maximum Likelihood Functions

optimization

r