Calcul de l'AIC «à la main» dans R

15

J'ai essayé de calculer l'AIC d'une régression linéaire dans R mais sans utiliser la AICfonction, comme ceci:

lm_mtcars <- lm(mpg ~ drat, mtcars)

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

Cependant, AICdonne une valeur différente:

AIC(lm_mtcars)
[1] 190.7999

Quelqu'un pourrait-il me dire ce que je fais mal?

luciano
la source
5
(sans vérifier encore votre réponse): Vous ne faites pas nécessairement quelque chose de mal, car la probabilité n'est en fait définie que jusqu'à une constante multiplicative; deux personnes peuvent calculer la log-vraisemblance et obtenir des nombres différents (mais les différences de log-vraisemblance sont les mêmes).
Glen_b -Reinstate Monica
1
La réponse de Hong Oois est liée à cette question, je pense. La formule utilisée par la fonction AICest -2*as.numeric(logLik(lm_mtcars))+2*(length(lm_mtcars$coefficients)+1).
COOLSerdash
luciano: Le "+1" dans cette formule à laquelle @COOLSerdash pointe provient du terme du paramètre de variance. Notez également que la fonction logLikindique que pour les lmmodèles, elle inclut «toutes les constantes» ... donc il y en aura log(2*pi)quelque part quelque part
Glen_b -Reinstate Monica
1
@Glen_b: Pourquoi dire que la vraisemblance n'est définie que jusqu'à une constante multiplicative? Après tout, lorsque vous comparez des modèles non imbriqués de différentes familles de distribution (par exemple avec AIC ou avec le test de Cox), vous devez vous souvenir de cette constante.
Scortchi - Réintégrer Monica
@Scortchi la définition n'est pas la mienne! Vous devrez le reprendre avec RAFisher. C'est ainsi depuis le début, je pense (1921). Qu'il est toujours défini de cette façon, au moins dans le cas continu, voir ici , par exemple, à la phrase commençant par «Plus précisément».
Glen_b -Reinstate Monica

Réponses:

18

Notez que l'aide sur la fonction logLikdans R indique que pour les lmmodèles, elle inclut «toutes les constantes» ... il y aura donc unlog(2*pi) in quelque part, ainsi qu'un autre terme constant pour l'exposant dans la vraisemblance. De plus, vous ne pouvez pas oublier de compter le fait que est un paramètre.σ2

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

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

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

AIC=2p2logL

mais notez que pour un modèle avec 1 variable indépendante, p = 3 (le coefficient x, la constante et )σ2

Ce qui signifie que c'est ainsi que vous obtenez leur réponse:

nrow(mtcars)*(log(2*pi)+1+log((sum(lm_mtcars$residuals^2)/nrow(mtcars))))
       +((length(lm_mtcars$coefficients)+1)*2)
Glen_b -Reinstate Monica
la source
Pourquoi dans votre calcul de divisez-vous uniquement par n et non par n - p ? s2nnp
Luke Thorburn du
1
Voir la définition de l' AIC: où le vecteur de paramètres, θ sont évalués au maximum ( à savoir tous les éléments de θ sont MLE); voir par exemple le critère d'information Wikipedia Akaike: Définition . Si vous n'êtes pas en divisant par n il y a dans le calcul de σ 2 , vous n'êtes pas le calcul de la MLE de σ 22logL(θ^)+2pθθ^nσ^2σ2 et donc pas vraiment calcul AIC - effet en vous seriez adaptez deux fois l'effet des paramètres d' ajustement. (Oui, beaucoup de gens le font mal)
Glen_b -Reinstate Monica
2logL=nlog(2π)+nlogsn+i(ei2/sn2)2πsn2
10

AIC2k2logLLknlogSrn+2(k1)Srn

Scortchi - Réintégrer Monica
la source