**From:** Achim Zeileis (*zeileis@ci.tuwien.ac.at*)

**Date:** Thu 26 Jun 2003 - 04:30:39 EST

**Next message:**Rafael Bertola: "[R] robust regression"**Previous message:**Edward Dick: "[R] logLik.lm()"**In reply to:**Edward Dick: "[R] logLik.lm()"**Next in thread:**Prof Brian Ripley: "Re: [R] logLik.lm()"

Message-id: <200306251830.h5PIUePF010374@thorin.ci.tuwien.ac.at>

On Wednesday 25 June 2003 20:23, Edward Dick wrote:

*> Hello,
*

*>
*

*> I'm trying to use AIC to choose between 2 models with
*

*> positive, continuous response variables and different error
*

*> distributions (specifically a Gamma GLM with log link and a
*

*> normal linear model for log(y)). I understand that in some
*

*> cases it may not be possible (or necessary) to discriminate
*

*> between these two distributions. However, for the normal
*

*> linear model I noticed a discrepancy between the output of
*

*> the AIC() function and my calculations done "by hand."
*

*> This is due to the output from the function logLik.lm(),
*

*> which does not match my results using the dnorm() function
*

*> (see simple regression example below).
*

*>
*

*> x <- c(43.22,41.11,76.97,77.67,124.77,110.71,144.46,188.05,171.18,
*

*>
*

*> 204.92,221.09,178.21,224.61,286.47,249.92,313.19,332.17,374.35) y <-
*

*> c(5.18,12.47,15.65,23.42,27.07,34.84,31.03,30.87,40.07,57.36,
*

*> 47.68,43.40,51.81,55.77,62.59,66.56,74.65,73.54)
*

*> test.lm <- lm(y~x)
*

*> y.hat <- fitted(test.lm)
*

*> sigma <- summary(test.lm)$sigma
*

*> logLik(test.lm)
*

*> # `log Lik.' -57.20699 (df=3)
*

*> sum(dnorm(y, y.hat, sigma, log=T))
*

*> # [1] -57.26704
*

*>
*

*> The difference in this simple example is slight, but
*

*> it is magnified when using my data.
*

That is because you are not using the ML estimate of the variance.

sigmaML <- sqrt(mean(residuals(test.lm)^2))

sum(dnorm(y, y.hat, sigmaML, log=T))

# [1] -57.20699

Best,

Z

*> My understanding is
*

*> that it is necessary to use the complete likelihood functions
*

*> for both the Gamma model and the 'lognormal' model (no constants
*

*> removed, etc.) in order to accurately compare using AIC.
*

*>
*

*> Can someone point out my error, or explain this discrepancy?
*

*>
*

*> Thanks in advance!
*

*>
*

*> ______________________________________________
*

*> R-help@stat.math.ethz.ch mailing list
*

*> https://www.stat.math.ethz.ch/mailman/listinfo/r-help
*

______________________________________________

R-help@stat.math.ethz.ch mailing list

https://www.stat.math.ethz.ch/mailman/listinfo/r-help

**Next message:**Rafael Bertola: "[R] robust regression"**Previous message:**Edward Dick: "[R] logLik.lm()"**In reply to:**Edward Dick: "[R] logLik.lm()"**Next in thread:**Prof Brian Ripley: "Re: [R] logLik.lm()"

*
This archive was generated by hypermail 2.1.3
: Tue 01 Jul 2003 - 09:12:09 EST
*