Как я могу получить стандартные ошибки для моих 4 параметров, когда матрица Гессе из solnp равна 5 на 5?

Я использую функцию solnp() в пакете R Rsolnp для решения нелинейной регрессии с ограничениями. Работает хорошо, сходится без проблем. Я хочу использовать матрицу Гессе для вычисления стандартных ошибок оценок четырех параметров, но гессиан не 4 на 4, как я ожидал, а 5 на 5. Я огляделся на SO и не увидел никого с неожиданным Гессенский размер. Все примеры, которые я нашел с напечатанными гессианами, показали, что они имеют ожидаемый размер p на p (например, 2x2, 3x3 и 4x4).

Как я могу получить стандартные ошибки для моих 4 параметров из этого гессиана 5 на 5?

df <- data.frame(
  Recruit.N = c(78.4, 79.8, 106, 57.4, 81.7, 94.4, 74.1, 42, 61.6, 47.7, 61.8, 
    28.1, 32.3, 19, 23.4, 20.1, 27), 
  Stock.5 = c(66.6, 90.3, 138.5, 79.8, 77.3, 78.4, 79.8, 106, 57.4, 81.7, 94.4, 
    74.1, 42, 61.6, 47.7, 61.8, 28.1), 
  Stock.6 = c(25.2, 66.6, 90.3, 138.5, 79.8, 77.3, 78.4, 79.8, 106, 57.4, 81.7, 
    94.4, 74.1, 42, 61.6, 47.7, 61.8), 
  Stock.7 = c(23.8, 25.2, 66.6, 90.3, 138.5, 79.8, 77.3, 78.4, 79.8, 106, 57.4, 
    81.7, 94.4, 74.1, 42, 61.6, 47.7)
)

lossfcn <- function(parz, mydat) {
  alpha <- parz[[1]]
  beta <- parz[[2]]
  p5 <- parz[[3]]
  p6 <- parz[[4]]
  p7 <- 1 - p5 - p6
  S <- with(mydat, p5*Stock.5 + p6*Stock.6 + p7*Stock.7)
  Obs <- mydat$Recruit.N
  Pred <- alpha * S * exp(-beta*S)
  Resid <- log(Obs) - log(Pred)
  sigma <- sqrt(mean(Resid^2))
  LL <- dlnorm(Obs, meanlog=log(Pred), sdlog=sigma, log=TRUE)
  -sum(LL)
}

inequal <- function(parz, mydat) {
  parz[3] + parz[4]
}

library(Rsolnp)
solnp(pars=c(1, 0.008, 1/3, 1/3), fun=lossfcn, mydat=df,
  ineqfun=inequal, ineqLB=0, ineqUB=1, 
  LB=c(0, 0, 0, 0), UB=c(1000, 1000, 1, 1), control=list(trace=0))

$pars
[1] 6.731317e-01 1.888572e-10 8.141363e-01 1.858631e-01

$convergence
[1] 0

$values
[1] 79.87150 75.50927 75.50927 75.50927

$lagrange
          [,1]
[1,] -2.028222

$hessian
           [,1]          [,2]         [,3]        [,4]         [,5]
[1,]  0.3350868 -3.359077e-01     17.84919  -0.4306057   -0.3382811
[2,] -0.3359077  1.993956e+02 -10161.63351  -7.0844295   -2.2749785
[3,] 17.8491854 -1.016163e+04 548099.69224 -85.9544831 -224.0362766
[4,] -0.4306057 -7.084429e+00    -85.95448  25.1086694    5.8817704
[5,] -0.3382811 -2.274979e+00   -224.03628   5.8817704    4.1978178

$ineqx0
[1] 0.9999995

$nfuneval
[1] 142

$outer.iter
[1] 3

$elapsed
Time difference of 0.03016496 secs

$vscale
[1] 1 1 1 1 1 1

person Jean V. Adams    schedule 23.01.2019    source источник


Ответы (1)


В отличие от трех сообщений, на которые вы ссылаетесь, у вас есть ограничение неравенства. Проверьте ineqx0 в возвращаемых значениях: в других сообщениях есть NULL, а у вас есть 0.9999995. С ограничением неравенства возникает незначительная переменная, поэтому проблема дополняется. Возвращаемая матрица Гессе предназначена для этого расширенного набора параметров. Просто сохраните первую подматрицу 4 x 4 hessian для нужных параметров.

person Zheyuan Li    schedule 26.01.2019
comment
Чтобы быть уверенным, что я выбрал правильную подматрицу, я изменил порядок параметров и обнаружил, что первая строка и столбец матрицы Гессе соответствуют дополнительному параметру неравенства. Итак, подматрица, соответствующая моим исходным четырем параметрам, равна hessian[-1, -1]. - person Jean V. Adams; 30.01.2019