Cálculo de AIC "a mano" en R

15

He intentado calcular el AIC de una regresión lineal en R pero sin usar la AICfunción, así:

lm_mtcars <- lm(mpg ~ drat, mtcars)

nrow(mtcars)*(log((sum(lm_mtcars$residuals^2)/nrow(mtcars))))+(length(lm_mtcars$coefficients)*2)
[1] 97.98786

Sin embargo, AICle da un valor diferente:

AIC(lm_mtcars)
[1] 190.7999

¿Alguien podría decirme qué estoy haciendo mal?

luciano
fuente
55
(sin verificar su respuesta todavía): no necesariamente está haciendo nada malo, ya que la probabilidad en realidad solo se define hasta una constante multiplicativa; dos personas pueden calcular la probabilidad logarítmica y obtener números diferentes (pero las diferencias en la probabilidad logarítmica son las mismas).
Glen_b -Reinstale a Monica
1
La respuesta de Hong Oois está relacionada con esta pregunta, creo. La fórmula que AICusa la función es -2*as.numeric(logLik(lm_mtcars))+2*(length(lm_mtcars$coefficients)+1).
COOLSerdash
luciano: El "+1" en esa fórmula @COOLSerdash apunta a que surge del término del parámetro de varianza. Tenga en cuenta también que la función logLikdice que para los lmmodelos incluye 'todas las constantes' ... por lo que habrá un log(2*pi)lugar allí en alguna parte
Glen_b -Reinstate Monica
1
@Glen_b: ¿Por qué decir que la probabilidad se define solo hasta una constante multiplicativa? Después de todo, al comparar modelos no anidados de diferentes familias de distribución (por ejemplo, con AIC o con la prueba de Cox), debe recordar esa constante.
Scortchi - Restablece a Monica
@Scortchi, ¡la definición no es mía! Tendrás que retomarlo con RAFisher. Ha sido así desde el principio, creo (1921). Que todavía se define de esa manera, al menos en el caso continuo, vea aquí , por ejemplo, en la oración que comienza 'Más precisamente'.
Glen_b -Reinstale a Monica

Respuestas:

18

Tenga en cuenta que la ayuda sobre la función logLiken R dice que para los lmmodelos incluye 'todas las constantes' ... por lo que habrá un log(2*pi)allí en algún lugar, así como otro término constante para el exponente en la probabilidad. Además, no puede olvidar contar el hecho de que es un parámetro.σ2

L(μ^,σ^)=(12πsn2)nexp(12i(ei2/sn2))

2logL=nlog(2π)+nlogsn2+i(ei2/sn2)

=n[log(2π)+logsn2+1]

AIC=2p2logL

pero tenga en cuenta que para un modelo con 1 variable independiente, p = 3 (el coeficiente x, la constante y )σ2

Lo que significa que así es como obtienes su respuesta:

nrow(mtcars)*(log(2*pi)+1+log((sum(lm_mtcars$residuals^2)/nrow(mtcars))))
       +((length(lm_mtcars$coefficients)+1)*2)
Glen_b -Reinstate a Monica
fuente
¿Por qué en tu cálculo de solo estás dividiendo entre n y no n - p ? s2nnp
Luke Thorburn
1
Véase la definición de AIC: donde el vector de parámetros, θ se evalúan en el máximo (es decir, todos los elementos de θ son MLEs); por ejemplo, ver el criterio de información de Wikipedia Akaike: Definición . Si usted no está dividiendo por n existe en el cálculo de σ 2 , usted no está calculando el MLE de σ 2 y así que no es el cálculo de la AIC - en efecto, que estaría ajustando dos veces por el efecto de los parámetros de ajuste. (Sí, mucha gente lo hace mal)2logL(θ^)+2pθθ^nσ^2σ2
Glen_b: reinstala a Monica el
2logL=nlog(2π)+nlogsn+i(ei2/sn2)2πsn2
10

La AICfunción da2k-2Iniciar sesiónL, dónde L es la probabilidad y k is the number of estimated parameters (including the intercept, & the variance). You're using nlogSrn+2(k1), where Sr is the residual sum of squares, & n is the sample size. These formulæ differ by an additive constant; so long as you're using the same formula & looking at differences in AIC between different models where the constants cancel, it doesn't matter.

Scortchi - Reinstate Monica
fuente