Maximum Likelihood Estimation

• How do we actually estimate parameters for generalized linear models?
  • Take our titanic example. You might think that we could simply take the logit transformations of our $y$'s and then use GLS to correct for heteroskedastacity.
  • The problem here is that our $y$'s only take the values of one and zero which will both be undefined in the logit transformation.
  • Therefore, we need a new technique - this is where maximum likelihood estimation comes in.

• Crash couse in Calculus
  • It helps to know a few basic things from calculus.
  • The derivative
    • Start with a function $f(x)$. This is some mathematical equation involving x.
      $$f(x)=5+28x-7x^{2}$$

      

    • Every function has a derivative. The derivative is basically the exact slope of the function at $x$.
      $$\frac{df(x)}{dx}=28-14x$$

    • Local minimums and maximums of the function $f(x)$ will occur when the slope is zero. (show graphically).
      

      $$= 28-14x$$ 0

      $$14x = 28$$

      $$x = 2$$

      
      
    • You can tell if it is a local minimum or maximum by looking at the second derivative, which is the derivative of the derivative. If this is positive at that point, then it is a minimum. If it is negative, then it is a maximum.
      $$\frac{d^{2}f(x)}{dx}=-14$$

    • So this is a maximum.

• The logic of maximum likelihood estimation is as follows:
  1. We have some process governed by some parameter(s), (say $\theta$) which generates a set of observed data ( $x_{1},x_{2},\ldots,x_{n}$ or something).
  2. We ask what is the likelihood, given the process, that we observe the given data? This leads to a likelihood function, $L(\theta)$.
  3. The trick is we don't know $\theta$. Rather it is usally the data we have and $\theta$ that we want to estimate. We choose a $\hat{\theta}$ such that it maximizes our likelihood function given the data.

• Example with a binomial proportion.
  • We have observed the outcome from a series of $n$ coin tosses. From this coin toss, we have observed $x$ heads. $n$ and $x$ are our observed data.
  • We know that the process of generating heads on coin tosses is a binomial process, because we only have two possible outcomes, and the probability $p$ is the same on each coin toss. $p$ (the probability of a heads on a single coin toss) is our unknown parameter.
  • Since this is a binomial process our likelihood of observing a particular number of heads in $n$ trials is given by the following likelihood function:
    $$L(p)=\binom{n}{x}p^{x}(1-p)^{n-x}$$
    (review this formula a bit)
  • We don't know the value of $p$, but we want to estimate it by choosing $\hat{p}$ that maximizes $L(p)$ for the given $n$ and $x$.
  • For particular values of $n$ and $p$, we could get $\hat{p}$ by trial and error, but that might take a long time. A computer might simplify it.

    Let's take the example, where $x=13$ and $n=50$. We can use the computer to calculate the likelihood functions for every $p$ at intervals of 0.01 from 0 to 1.

    

    The $p$ which would maximize our likelihood of actually observing 13 heads on 50 tosses is about 0.26 - this is what we would expect by common sense.

  • Even so, we would like to get a general result that can be expressed in the abstract terms $n$ and $p$.
  • Thus, we need to find the maximum of the likelihood function for $p$. We can do this using calculus. The point where the derivative of the likelihood function is zero and the second derivative is negative is the maximum.
  • It turns out, that working with the log of the likelihood function is usually easier than working with the likelihood function itself. The results for the log-likelihood function hold for the likelihood function (the maximum is the same), so logging it does no damage.
    $$\log L(p)=\log \binom{n}{x}+x\log(p)+(n-x)\log(1-p)$$

  • We take the derivative of this function with respect to $p$.
    $$\frac{\partial \log L(p)}{\partial p}=\frac{x}{p}-\frac{n-x}{1-p}$$

  • When this function equals zero, we will have either a minimum or a maximum. So solve for $p$.
    

    $$= \frac{x}{p}-\frac{n-x}{1-p}$$ 0

    $$\frac{n-x}{1-p} = \frac{x}{p}$$

    $$p(n-x) = x(1-p)$$

    $$pn-px = x-px$$

    $$pn = x$$

    $$p = \frac{x}{n}$$

    
    
  • It turns out the second derivative at this point is negative (not shown), so this is indeed a maximum. Thus our best estimate of $p$ is $\hat{p}=x/n$. The same as our intuition would have told us.

• This example was relatively straightforward. But we can use the same technique with the linear model.
• Our parameters are the $\beta$'s, which relate our $x$ variables to our $y$ variable.
• For some sample, we have observed data on the $y$ variable and the set of $x$ variables. Given this data, there are some $\beta$ values which would have maximized our likelihood of observing it. By MLE, these values are our estimates of the $\beta$ parameters.
• How do we figure out what $\beta$'s would have maximized our likelihood of observing the actual data? That is complicated. The remainder of appendix B provides a good introduction to the general technique, but it is highly technical. All we need to know for now is that computers are very good at finding such things. We just need to know intuitively what process guides their decision.

• The use of MLE gives us a new method of model selection.
  • We can think of the actual log-likelihood returned from maximizing the likelihood function as a measure of how well it fits.
  • However, we can't use the log-likelihood $Log L$ directly, because different sample sizes will have different likelihoods.
  • What we can use is what is called deviance or $G^{2}$, which measures how much our current model deviates from the saturated model. If $L_{f}$ is the likelihood of the saturated model, then:
    $$G^{2}=-2 \log (L_{c})/L_{f})=-2(\log L_{c}-\log L_{f})$$
    The log of the likelihood will always be negative or zero, so multiplying it by a negative number will create a positive number. The closer this deviance is to zero, the closer our likelihood is to the saturated model.
  • For most of the individual-level data we work with the likelihood of the saturated model is 1, so $\log L_{f}=0$, so our deviance simplifies:
    $$G^{2}=-2 \log L_{c}$$

  • We can use deviances to compare two nested models using the Likelihood Ratio Test (LRT).
  • Differens between nested models follow our friend the $\chi^{2}$ distribution. Therefore, if we have two models $n$ and $u$ in which $n$ is nested within $u$, the value $G^{2}_{u}-G^{2}_{n}$ has a $\chi^{2}$ distribution with DF of $df_{u}-df_{n}$, assuming model $u$ is the correct model. This method can be used to compare two or more nested models.
  • You can also use deviances to generate $BIC$ statistics. To generate a BIC against the saturated model:
    $$BIC=G^{2}-df\log n$$
    In most cases, however we will want BIC' relative to the null model. To calculate this first get the deviance for the null model (often reported in statistical software). Then:
    $$BIC'=(G^{2}_{0}-G^{2})+p\log n$$
    where p is the number of parameters fit in the current model.
  • Although BIC and BIC' have different interpretations directly, when you compare the BIC of BIC' of two models, you should get the same difference.