Namans
Namans

Reputation: 35

How to find the maximum likelihood estimate of p in a binomial distribution characterized by 9 successes in 20 trials using R?

Find the maximum likelihood estimate of p in a binomial distribution characterized by 9 successes in 20 trials. Show graphically that this is the maximum. Include your R code with your answers.

This is a prompt I've been given for a homework assignment but the teacher never explained how to do it. I understand the basic concept, but I'm not sure how to find the estimate of p in R, or how to graph it. Any advice/help would be greatly appreciated!

Upvotes: 0

Views: 76

Answers (1)

Rui Barradas
Rui Barradas

Reputation: 76402

The exercise is not that difficult.

  1. Write down two functions, the textbook likelihood and log-likelihood;
  2. Use optim to maximize them;
  3. Below the tol argument is smaller than the default. Since both functions are continuous, finite and convex their maxima are guaranteed to exist and a smaller tol will find better maxima values;
  4. Finally, plot the functions and the maxima found earlier.
ll <- function(x, n, k) choose(n, k) * x^k * (1 - x)^(n - k)

log_ll <- function(x, n, k) lchoose(n, k) + k*log(x) + (n - k)*log(1 - x)

xmax <- optimize(ll, c(0, 1), n = 20, k = 9, maximum = TRUE, tol = .Machine$double.eps^0.5)
xmax$maximum
#> [1] 0.45

xmax_log <- optimize(log_ll, c(0, 1), n = 20, k = 9, maximum = TRUE, tol = .Machine$double.eps^0.5)
xmax_log$maximum
#> [1] 0.45

# save the default graphics parameters
old_par <- par(mfrow = c(2, 1), mai = c(0.8, 1, 0.1, 1))
#
curve(ll(x, n = 20, k = 9), from = 0, to = 1, xlab = "")
segments(x0 = xmax$maximum, y0 = -1, y1 = xmax$objective, lty = "dashed")
segments(x0 = -1, x1 = xmax$maximum, y0 = xmax$objective, lty = "dashed")
points(xmax$maximum, xmax$objective, col = "red", pch = 16)
#
curve(log_ll(x, n = 20, k = 9), from = 0, to = 1, xlab = "Binomial proportion")
segments(x0 = xmax_log$maximum, y0 = -50, y1 = xmax_log$objective, lty = "dashed")
segments(x0 = -1, x1 = xmax_log$maximum, y0 = xmax_log$objective, lty = "dashed")
points(xmax_log$maximum, xmax_log$objective, col = "red", pch = 16)

#
# restore the default graphics parameters
par(old_par)

Created on 2022-10-11 with reprex v2.0.2

Upvotes: 1

Related Questions