Let's say I have an outcome that is exponentially distributed, so $p(y|\lambda_i) = \lambda_i e^{-\lambda_i y}$. However, we also know that $\log{\lambda_i} = \beta_0+\beta_1x_i$. If we want to estimate the MLE's for $\beta_0$ and $\beta_1$, what's the best way to do it given data?
I was thinking of using the following code:
summary(glm(y~x,family=Gamma(link="log")))
Not 100% sure what the various arguments are going though.
This has been answered on the R help list by Adelchi Azzalini: the important point is that the dispersion parameter (which is what distinguishes an exponential distribution from the more general Gamma distribution) does not affect the parameter estimates in a generalized linear model, only the standard errors of the parameters/confidence intervals/p-values etc.; in R an estimate of the dispersion parameter is automatically reported, but as Azzalini comments, summary.glm allows the user to specify the dispersion parameter. So, as stated by Azzalini,
The Gamma family is parametrised in glm() by two parameters: mean and dispersion; the "dispersion" regulates the shape.
So [one] must fit a GLM with the Gamma family, and then produce a "summary" with dispersion parameter set equal to 1, since this value corresponds to the exponential distribution in the Gamma family.
In practice:
fit <- glm(formula =..., family = Gamma(link="log"))
summary(fit,dispersion=1)
[Azzalini had family=Gamma, i.e. using the default inverse link; I changed it to specify the log link as in your question.]
Using a GLM call as you suggest there is the easiest correct approach, but to actually make the Gamma into an exponential you can specify the dispersion to be 1. It will not change the fitted mean / coefficients, but it impacts the standard errors.
[i.e. Your suggested call of summary(glm(y~x,family=Gamma(link="log"))) should give you what you want, but if you're interested in significance of coefficients and so on under the exponential assumption, you'd add ,dispersion=1 before the final parenthesis. If you want to fit Gamma GLMs more generally, there are several useful helper functions in the package MASS that comes with R but is not loaded by default.]
(Another alternative might be to consider using parametric survival models which also offer ways to fit exponential, Weibull and various other models.)
bbmle::mle2will do this. You could try to derive and solve the score equations but I don't know how far you'd get ... $\endgroup$ – Ben Bolker Oct 16 '16 at 0:34bbmle::mle2(y~dexp(exp(lograte)), parameters=list(lograte~x), start=list(logmu=0), data=...)should do it, I think ... $\endgroup$ – Ben Bolker Oct 16 '16 at 0:46