r - 对于带有bobyqa方法的optimx()，如何设置初始值

Question

我试图用 R 包“optimx”来最大化我的可能性。这是我的代码。使用初始值 (5,5) 和 (1,1)，我得到了不同的最大化似然。我也尝试过不同的方法，如“Nelder=Mead”，但不同方法下估计的对数似然度不同......

library('optimx')
n=225
X = matrix(runif(225),ncol=1)
e2 = matrix(runif(225,0,2),ncol=1)

set.seed(123)

这是生成我将使用的一些数据的函数

get_mls_basis<- function(p){

  depth <- ceiling(runif(1)*p)

  knot <- matrix(rep(0,depth+1),ncol=1)

  lr <- runif(1) > 0.5

  x <- matrix(rep(0,n),ncol=1)

  not_finished <- 1

  while (not_finished == 1) {
    data_indx = ceiling(runif(1)*n)
    var = matrix(rep(0,depth),ncol=1) 
    for (j in 1:depth) {
      not_ok <- 1
      while (not_ok == 1) {
        ind <- ceiling(runif(1)*p)
        if (!is.element (ind,var[1:j]))
        {
          var[j] <- ind
          not_ok <- 0
        }
      }
    }

    x_v <- as.matrix(X[data_indx, var])
    knot[1:depth] <- rgamma(depth,1,1)
    knot[1:depth] <- knot[1:depth] / sqrt(sum(knot^2))
    knot[depth+1] <- -x_v %*% knot[1:depth]

    ones <- matrix(rep(1,n),ncol=1)

    temp <- as.matrix(cbind(X[,var], ones)) %*% knot

    if (lr == 0) {
      for (i in 1:n)
      {
        temp[i] <- max(0,temp[i])
      }
    }
    else {
      for (i in 1:n)
      {
        temp[i] <- min(0,temp[i])
      }
    }

    x <- temp
    not_finished <- all(x==0)
  }

  mx <- mean(x)
  stx <- sd(x)
  x <- (x-mx)/stx
  x



}

这是我的日志可能性

Lik1<-function(theta, basis){
  theta0=theta[1]
  theta1=theta[2]
  L=-n/2*log(theta0)-sum(basis/2)*log(theta1)-0.5/theta0*sum(e2/theta1^basis)
  return(L)
}
 basis1=get_mls_basis(1)

这里我使用 5 作为初始值

  optimx(par=c(5,5), Lik1,
         basis=basis1,method='bobyqa',control = list(maximize=TRUE))

标签： rmle