I'm fitting models to data with the
mle2() function in the bbmle package.
In an example problem (reproducible code below) I am fitting a count response variable to a single predictor with a linear relationship (polar bears, PB, and fraction of coastline filled with sea ice, Ice). I created a simulated dataset by sampling from a Poisson distribution for polar bears using a linear function of sea ice and an intercept and then fit the model with
mle2(). Since this is a linear model I can fit the same model with glm(), of course. The two approaches (glm and mle2) provide identical results only if one uses the identity link function, rather than the (canonical) log link (or any other link).
My question is: should I be using a log link in the mle2() fitting (and appropriate link functions for other distributions) when using mle2()?
I understand that fitting non-Gaussian glm() models often uses iterated weighted least squares, and that canonical link functions often translate the expected value of the linear model onto the -infinity to +infinity range (e.g. a log link takes a distribution that ranges from 0 to +infinity and translates it to -infinity to +infinity), and so I can see how a link function would aid in the fitting process. However, I wouldn’t have thought there was anything wrong with the mle2() fit, which is equivalent to using the identity link. But given that using glm() with the canonical link and mle2() produce different answers I'm uncertain about which one is preferred and why. Thank you. marm
require(bbmle) set.seed(67) Ice=runif(20,0,1); PB=rpois(20,lambda=20+50*Ice);plot(PB~Ice) d=data.frame(PB=PB,Ice=Ice) F1=mle2(PB~dpois(lambda=c0+c1*Ice),start=list(c0=10,c1=30),data=d);summary(F1) F3=glm(PB~Ice,family=poisson(link=identity),data=d);summary(F3) #gives identical results to F1 F4=glm(PB~Ice,family=poisson,data=d);summary(F4) #canonical link - gives different results for predicted values (and obviously for coefs) points(predict(F3)~Ice,col="blue") #predicted values from identity link or mle2() fit points(predict(F4,type="response")~Ice,col="red") #predicted values from canonical link - curvature due to log link is apparent in differences from F3 predictions