7 votos

¿Por qué la fórmula del AIC en R parece utilizar un parámetro más de lo esperado?

Utilizaré un ejemplo para que puedas reproducir los resultados

# mortality 
mort = ts(scan("http://www.stat.pitt.edu/stoffer/tsa2/data/cmort.dat"),start=1970, frequency=52)

# temperature
temp = ts(scan("http://www.stat.pitt.edu/stoffer/tsa2/data/temp.dat"), start=1970, frequency=52)

#pollutant particulates
part = ts(scan("http://www.stat.pitt.edu/stoffer/tsa2/data/part.dat"), start=1970, frequency=52)

temp = temp-mean(temp)
temp2 = temp^2
trend = time(mort)

Ahora, ajuste un modelo para los datos de mortalidad

fit = lm(mort ~ trend + temp + temp2 + part, na.action=NULL)

Lo que quiero ahora es reproducir el resultado del comando AIC

AIC(fit)
[1] 3332.282

Según el archivo de ayuda de R para AIC, AIC = -2 * log.likelihood + 2 * npar. Si estoy en lo cierto, creo que log.likelihood se da usando la siguiente fórmula:

n = length(mort)
RSS = anova(fit)[length(anova(fit)[,2]),2] # there must be better ways to get this, anyway
(log.likelihood <- -n/2*(log(2*pi)+log(RSS/n)+1))

 [1] -1660.135

Esto es aproximadamente igual a

logLik(fit)
'log Lik.' -1660.141 (df=6)

Por lo que veo, el número de parámetros del modelo es de 5 (¿cómo puedo obtener este número programáticamente?). Así que el AIC debe ser dado por:

-2 * log.likelihood + 2 * 5
[1] 3330.271

Ooops, parece que debería haber utilizado 6 en lugar de 5 como número de parámetros. ¿Qué hay de malo en esos cálculos?

13voto

Patrick Puntos 183
> -2*logLik(fit)+2*(length(fit$coef)+1)
[1] 3332.282

(lo has olvidado; tienes 6 parámetros porque $\sigma_{\epsilon}$ ¡también tiene que ser estimado!

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X