Em R, dado um resultado de optim com uma matriz de hessiana, como calcular intervalos de confiança de parâmetro usando a matriz de hessiana?

22

Dado um resultado do otim com uma matriz hessiana, como calcular os intervalos de confiança dos parâmetros usando a matriz hessiana?

fit<-optim(..., hessian=T)

hessian<-fit$hessian

Estou interessado principalmente no contexto da análise de máxima verossimilhança, mas estou curioso para saber se o método pode ser expandido além.

Etienne Low-Décarie
fonte
2
Esta pergunta é muito vaga. Que tipo de intervalos de confiança? Em que tipo de modelos você está interessado? Considere dar uma olhada na literatura antes de postar 3 perguntas seguidas.
O método deve ser independente do intervalo e modelo.
Etienne Low-Décarie
Que função você otimiza?
Stéphane Laurent
Foi-me dito que eu deveria ser capaz de fazer isso independentemente do modelo que eu uso.
Etienne Low-Décarie
Corey Chivers forneceu a resposta.
Etienne Low-Décarie

Respostas:

32

Se você está maximizando uma probabilidade, a matriz de covariância das estimativas é (assintoticamente) o inverso do negativo do hessiano. Os erros padrão são as raízes quadradas dos elementos diagonais da covariância ( de outros lugares da Web ! , do Prof. Thomas Lumley e Spencer Graves, Eng.).

Para um intervalo de confiança de 95%

fit<-optim(pars,li_func,control=list("fnscale"=-1),hessian=TRUE,...)
fisher_info<-solve(-fit$hessian)
prop_sigma<-sqrt(diag(fisher_info))
prop_sigma<-diag(prop_sigma)
upper<-fit$par+1.96*prop_sigma
lower<-fit$par-1.96*prop_sigma
interval<-data.frame(value=fit$par, upper=upper, lower=lower)

Observe que:

  • Se você está maximizando o log (probabilidade), o ponto negativo do hessian é a "informação observada" (como aqui).
  • Se você MINIMIZAR um "desvio" = (-2) * log (probabilidade), então a METADE do hessiano é a informação observada.
  • No caso improvável de maximizar a probabilidade em si, é necessário dividir o negativo do hessiano pela probabilidade de obter as informações observadas.

Veja isso para obter outras limitações devido à rotina de otimização usada.

Etienne Low-Décarie
fonte
Corey Chivers forneceu a resposta.
Etienne Low-Décarie
2
(+1) O inverso do hessiano negativo é um estimador da matriz de covariância assintótica - eu sei que isso aparece no seu código, mas acho importante ressaltar.
Macro
1
Excelente resposta, os limites superior e inferior devem ser lidos upper<-fit$par+1.96*(prop_sigma/sqrt(n)) lower<-fit$par-1.96*(prop_sigma/sqrt(n))? Obrigado
forecaster
3
Por que não excluir a linha 4?
Jason
2
A inclusão da linha é prop_sigma<-diag(prop_sigma)um erro?
Mark Miller