# Let $X_1,\dots, X_n$ be random sample from $Bernoulli(p)$. Which estimator is better?

Let $$X_1,\dots, X_n$$ be random sample from $$Bernoulli(p)$$. Compare the risks of the squared loss of two estimators of $$p$$:

$$\hat{p}_1=\bar{X}, \, \hat{p}_2=\frac{n\bar{X}+\alpha}{\alpha+\beta+n}$$ where $$\bar{X}$$ is the sample mean, and $$\alpha, \beta$$ are two positive constants. Which estimator is better?

I compute the MSE of these two estimators but I have no idea how to compare them.

I get $$MSE(\hat{p}_1)=Var(\hat{p}_1)+(Bias(\hat{p}_1))^2=\frac{p(1-p)}{n}$$ and $$MSE(\hat{p}_2)=Var(\hat{p}_2)+(Bias(\hat{p}_2))^2=\frac{np(1-p)}{(\alpha+\beta+n)^2}+\frac{(\alpha-p(\alpha+\beta))^2}{(\alpha+\beta+n)^2}$$

Then I compare the difference between them: $$MSE(\hat{p}_1)-MSE(\hat{p}_2)=\frac{p(1-p)(\alpha+\beta)(\alpha+\beta+2n)-n(\alpha-p(\alpha+\beta))^2}{(\alpha+\beta+n)^2}$$

• Idea: analyse the following inequality: $MSE(\hat{p}_1) \geq MSE(\hat{p}_2)$. Investigate for which $p$, $n$, $\alpha$ and $\beta$ the inequality holds. You might want to try and fix either $p$ or $n$ at first if you find working with both $p$ and $n$ as variables to be tricky Commented Feb 27, 2023 at 7:13

You have been given the MLE as well as the posterior mean under a Beta prior (see e.g. here).

As you write (in slightly different notation with $$\theta=p$$, $$k=n\bar X$$, $$\alpha_0=\alpha$$, $$\hat \theta_1=\hat p_1$$, $$R=MSE$$ etc.), \begin{align*} R(\theta, \hat{\theta}_2) =& \, V_\theta (\hat{\theta}_2) + ({bias}_\theta (\hat{\theta}_2))^2 \\ =& \, V_\theta \left( \frac{k + \alpha_0}{\alpha_0 + \beta_0 + n} \right) + \left( E_\theta \left( \frac{k + \alpha_0}{\alpha_0 + \beta_0 + n} \right) - \theta \right)^2 \\ =& \, \frac{n\theta(1 - \theta)}{(\alpha_0 + \beta_0 + n)^2} + \left( \frac{n\theta + \alpha_0}{\alpha_0 + \beta_0 + n} - \theta \right)^2. \end{align*} Which estimator then is "better" indeed depends on the prior parameters, sampls size as well as the true $$\theta$$. This is intuitive as a prior that puts lots of probability on regions of the parameter space where the true value happens to be can be expected to produce a good Bayesian estimator.

A well-known and interesting choice is to let $$\alpha_0 = \beta_0 = \sqrt{0.25n}$$. The resulting estimator is $$\hat{\theta}_2 = \, \frac{k + \sqrt{0.25n}}{n + \sqrt{n}}$$ Then, the risk function is independent of $$\theta$$: $$R(\theta,\hat{\theta}_2) = \, \frac{n}{4(n + \sqrt{n})^2},$$ and you can easily find the region in which $$\hat \theta_1$$ is outperformed in terms of maximum risk, and that that region shrinks with $$n$$ (that it shrinks with $$n$$ is intuitive, the sample average being the asymptotically efficient MLE).

The prior is "well-known" because there is a theorem (which I think I have come across in a book by Wasserman, but surely is also available elsewhere) saying that if $$\hat{\theta}$$ is a Bayes rule (a decision rule that minimizes the Bayes risk, e.g. the integral over $$\theta$$ for the MSE for squared loss with respect to some prior) with respect to some prior $$\pi$$ and if $$\hat{\theta}$$ has a constant risk \begin{align*} R(\theta, \hat{\theta}) = c \end{align*} for some $$c$$, then $$\hat{\theta}$$ is minimax.

A rule is minimax if, in a class of estimators $$\tilde{\theta}$$, it minimizes the maximum risk.

Formally, $$\hat{\theta}$$ is minimax if $$\sup\limits_{\theta} R(\theta, \hat{\theta}) = \, \inf\limits_{\tilde{\theta}} \sup\limits_{\theta} R(\theta, \tilde{\theta}).$$

In this example, \begin{align*} \sup\limits_{\theta} R(\theta, \hat{\theta}_1) =& \, \max\limits_{0 \leq \theta \leq 1} \frac{\theta(1-\theta)}{n} = \frac{1}{4n} \end{align*} and \begin{align*} \sup\limits_{\theta} R(\theta, \hat{\theta}_2) =& \, \max\limits_{\theta} \frac{n}{4(n + \sqrt{n})^2} = \frac{n}{4(n + \sqrt{n})^2} \end{align*} and we observe that $$\sup\limits_{\theta} R(\theta, \hat{\theta}_2)<\sup\limits_{\theta} R(\theta, \hat{\theta}_1)$$ On the other hand, consider a uniform prior $$\pi(\theta) = 1$$. Then, the Bayes risks are \begin{align*} r(\pi, \hat{\theta}_1) =& \, \int R(\theta, \hat{\theta}_1) \, d\theta = \int \frac{\theta(1-\theta)}{n} \, d\theta = \frac{1}{6n} \end{align*} and \begin{align*} r(\pi, \hat{\theta}_2) =& \, \int R(\theta, \hat{\theta}_2) \, d\theta = \frac{n}{4(n + \sqrt{n})^2}. \end{align*} For $$n \geq 20$$, $$r(\pi, \hat{\theta}_2) > r(\pi, \hat{\theta}_1)$$ (in the figure below, only the red dots are above the red dahes, with the ranking reversed for the other $$n$$) which suggests that $$\hat{\theta}_1$$ is a better estimator according to this criterior and this prior.

Schematically:

n <- c(15, 25, 30, 100)

theta <- seq(0.001, 0.999, 0.001)
risk.MLE <- sapply(n, function(n) theta*(1-theta)/n)
risk.posteriormean <- n/(4*(n + sqrt(n))^2)
Bayesrisk.MLE. <- 1/(6*n)

matplot(theta, risk.MLE, type="l", col=c("red", "blue", "green", "orange"), lty=1, lwd=2)
abline(h=risk.posteriormean, lwd=2, col=c("red", "blue", "green", "orange"), lty=2) # long dashes for the (Bayes) risk of the posterior mean
abline(h=Bayesrisk.MLE., lwd=2, col=c("red", "blue", "green", "orange"), lty=3) # dots for Bayes risk of MLE


So to make a long answer short: it depends!

• I like your answer (+1)! I am curious to know if there is a rationale or a reference for such a well-known choice. Commented Feb 27, 2023 at 8:48
• Fair question, I made an edit! Commented Feb 27, 2023 at 9:02
• I wonder whether the reverse of the theorem is also true. If $\hat\theta$ is minimax, then $R(\theta,\hat\theta)$ is constant. I can imagine a case where it is only constant in some part but not everywhere. Commented Feb 27, 2023 at 9:54
• My statistical decision theory knowledge is about to be exhausted (fwiw, the theorem as stated does not say the two properties are equivalent) - time for @Xi'an to weigh in :-)! Commented Feb 27, 2023 at 9:56
• If $\Lambda$ is a prior on $\Theta$ and $\int_\Theta R(\theta, \hat \theta) \mathrm d \Lambda(\theta)= \sup_\theta R(\theta, \hat \theta)$ holds (which is true for cst. risk) for a Bayes estimator $\hat \theta$ w.r.t. $\Lambda$, then $\hat \theta$ is minimax since then $\sup_\theta R(\theta,\tilde \theta) \geq \int_\Theta R(\theta,\tilde \theta) \mathrm d \Lambda(\theta) \geq \int_\Theta R(\theta,\hat \theta) \mathrm d \Lambda(\theta) = \sup_\theta R(\theta,\hat \theta)$. Clearly, if additionally $\hat \theta$ is unique Bayes, it is also unique minimax (replace the second $\geq$ with $>$). Commented Feb 27, 2023 at 12:16