# Problems to obtain the maximum likelihood estimator

I'm trying to calcul the MLE from a Beta-Binomial distribution. However, I'm having problems defining the estimator.

I'm using the following function:

$$\widehat{\ell\,}(\theta\,;x)=\sum_{i=1}^n \ln f(x_i\mid\theta),$$ I tried to calcul the function and then apply the log() function in each point:

n = 9
k = seq(0,n)

y = VGAM::dbetabinom(x = k, size = n, prob = 0.09894759, rho = 0.08603263)

mle = sum(log(y))



The problem is that the results I get are different from when compared to some function already existing in R.

I'm using the function below and trying to estimate alpha e beta from it. The parameters n and k are predefined.

$$f(k\mid n,\alpha,\beta) = \frac{\Gamma(n+1)}{\Gamma(k+1)\Gamma(n-k+1)} \frac{\Gamma(k+\alpha)\Gamma(n-k+\beta)}{\Gamma(n+\alpha+\beta)} \frac{\Gamma(\alpha+\beta)}{\Gamma(\alpha)\Gamma(\beta)}$$

I imagine that I have a problem with the concept, but I can not find the correct result. I don't undersand why when I sum all log() values I don't get the MLE.

How do I make the correct calculation of this estimator?

• Where do you optimise the (log-)likelihood? Jul 1, 2019 at 15:28
• Actually, I'm trying to optimize this function with Excel's Solver, but when I compare my results to those obtained by functions in R, they are different. Jul 2, 2019 at 7:32

You are quite confused. I see two major confusions: you are confusing the idea of the log-likelihood with that of the log-density and you do not understand the idea of an estimator.

The density function of a beta binomial is given by the formula you listed and can be evaluated in R by the function dbetabinom although note that dbetabinom actually uses a different parameterization than the formula you listed. You can get the same parameterization if you use dbetabinom.ab.

The log-likelihood is a function of the parameters ($$\alpha$$, and $$\beta$$) given some observed data. So suppose you have 1000 independent observations, $$k_1, \ldots,k_{1000}$$ and you assume each one came from a BetaBinomial with the same $$\alpha$$ and $$\beta$$. This is usually written $$k_i \sim \text{BetaBinomial}(n, \alpha, \beta)$$

In R we can draw such a sample of independent observations from the BetaBinomial like this (sorry I used y in my code rather than k):

library(VGAM)
n = 10
alpha <- 600
beta <- 400

y <- rbetabinom.ab(1000, size=n, alpha, beta)

hist(y)


Here is the histogram of those observations: The log-likelihood, $$\ell(\alpha, \beta)=\log{L(\alpha, \beta)}$$ is the log of the likelihood function which is a function of the parameters. Because we assumed each observation is independent we can factorize this likelihood as

$$L(\alpha, \beta) = \prod_{i=1}^{1000} p(k_i|n, \alpha, \beta)$$ where $$p(k_i|n,\alpha,\beta)$$ is the formula you listed.

The MLE (maximum likelihood estimator) is the $$\alpha$$ and $$\beta$$ which maximize this log-likelihood function

$$(\hat{\alpha}, \hat{\beta}) = \displaystyle\text{argmax}_{\alpha, \beta}{\ell(\alpha, \beta)}$$

Now for many distributions you can get a closed form solution for the MLE. There is no closed form solution for the MLE of a beta-binomial. So you would have to obtain the MLE using computational means (Newton-Raphson algorithm etc).

Now for illustrative purposes Lets look at the log-likelihood function with one of parameters held constant at its true value:

f1 <- function(a){
sapply(a, function(ax){ sum(log(dbetabinom.ab(y, size=n, ax, beta)))})
}

f2 <- function(b){
sapply(b, function(bx){
sum(log(dbetabinom.ab(y, size=n, alpha, bx)))
})
}

par(mfrow=c(1,2))
curve(f1(x), 0, 1000)
abline(v=alpha, col="gray", lty=2)
curve(f2(x), 0, 1000)
abline(v=beta, col="gray", lty=2)
par(mfrow=c(1,1)) f1 is $$\ell(\alpha, \beta=400)$$ and similarly f2 is $$\ell(\alpha=600, \beta)$$ I have plotted the log-likelihood as these parameters vary and in a vertical gray dotted line I have the true value. You can see this true value is close to the value which would maximize these functions.

• I understood your example. What I'm trying to do is optimize this function using Excel's Solver. For this I need to select the objective cell and the variable cells. I already did the calculation using the RSS estimator (I select the RSS value as the objective and the alpha and beta parameters as variables). However, I can not define a cell with the MLE value to use as an objective cell in Excel. Jul 2, 2019 at 7:45