3
$\begingroup$

Let us have a sample $X_1, X_2, ... , X_n$ with $B_{p,m}$ distribution.

How to make an estimator for $m$ using MLE ? For simplicity, let us $n=1$

I am a little bit stuck, because I have a derivative of binominal coefficient which is uncommon in usual undergrad textbooks.. So, probably, I am doing something wrong.

Update: $p$ is known

$\endgroup$
6
  • 1
    $\begingroup$ You cannot take a derivative since $m$ is an integer. $\endgroup$
    – Xi'an
    Feb 28, 2021 at 11:21
  • $\begingroup$ The case with $p$ unknown $\endgroup$ Feb 28, 2021 at 13:27
  • 1
    $\begingroup$ With a sample size $n=1,$ the MLE should be one of the last things you consider using, because its theory relies on asymptotic behavior (and $n=1$ is as far from asymptotic as one can get!). In this case it is useful to consider what your loss function might be and proceed from there to select an appropriate estimator. $\endgroup$
    – whuber
    Feb 28, 2021 at 17:57
  • $\begingroup$ Can you confirm you have a single sample with only one element $n=1$? $\endgroup$ Mar 1, 2021 at 19:42
  • $\begingroup$ If you have just a single element PMF is constant and equal to $p$ so you cannot derive anything. $\endgroup$ Mar 1, 2021 at 19:49

3 Answers 3

1
$\begingroup$

We have the likelihood of $m$ by $$f(m)=\frac{m!}{(m-X)!X!}p^X(1-p)^{m-X}$$

Here $m\geq X$.

Now, let's consider when $\frac{f(m+1)}{f(m)}$ is less than 1. When it starts to be less than 1, we know that function $f$ reached a (local) maximum. After calculation we find

$$r = \frac{f(m+1)}{f(m)} = \frac{(m+1)(1-p)}{m+1-X},$$ and when $n\geq\left\lfloor\frac{X}{p}\right\rfloor$, we have $r\leq 1$. Thus, we conclude that $$\left\lfloor\frac{X}{p}\right\rfloor\text{ is the MLE of }m.$$

Here, $\lfloor y\rfloor$ denotes the largest integer that is less than $y\in\mathbb{R}.$

----- update -----

By the way, BruceET's numerical example confirms the conclusion.

$\endgroup$
6
  • 1
    $\begingroup$ Typo: when $m \to X$ or $m \to \infty$ we have $f(m) \to 0$. Also, you answer assumes that $p$ is known. In the general case where both $n$ and $p$ are unknown, the MLE of $n$ may be $\infty$. This corresponds to the cas where the Poisson lilkelihood is greater than the binomial likelihood. See my answer and the reference therein. $\endgroup$
    – Yves
    Feb 28, 2021 at 8:20
  • 3
    $\begingroup$ When $m=X$, $f(X)=p^X$, which is neither $0$ nor $-\infty$. $\endgroup$
    – Xi'an
    Feb 28, 2021 at 11:27
  • $\begingroup$ Thanks, updated. $\endgroup$
    – Tan
    Feb 28, 2021 at 16:18
  • $\begingroup$ You are using $X$ index to denote the number of elements of a sample which is wrong since you already used the $n$ as index to denote the cardinality. $\endgroup$ Mar 1, 2021 at 19:47
  • $\begingroup$ @EasyPoints What is the "number of elements of a sample"? $X$ means the number of successes in $m$ independent Bernoulli trials. $\endgroup$
    – Tan
    Mar 1, 2021 at 20:11
2
$\begingroup$

This is along the lines you started to try.
$\sum{X_i}$ is binomial with parameters $p$ and $mn$.
The likelihood function is $$f(m)$$ $$=\frac{\Gamma(m+n+1)}{\Gamma(m+n-\sum{X_i}+1)\Gamma(\sum{X_i}+1)}p^{\sum{X_i}}(1-p)^{mn-\sum{X_i}}$$ for $\max{X_i}\le m$.

Start by treating the likelihood function as if it were a function of a continuous variable $m \in (\max{X_i},\infty)$.

The derivative of the log of the likelihood is $$f'(m)=\log{(1-p)}+\psi (mn+1)-\psi{\left(mn-\sum{X_i}+1\right)}$$ where $\psi$ is the digamma function.
This simplifies to $$f'(m)=\log{(1-p)}+\frac{1}{mn-\sum{X_i}+1}+...+\frac{1}{mn}$$ The second derivative is $$f''(m)=\psi^{(1)} (mn+1)-\psi^{(1)}{\left(mn-\sum{X_i}+1\right)}$$ where $\psi^{(1)}$ is the polygamma function.
This simplifies to $$f''(m)=-\frac{1}{\left(mn-\sum{X_i}+1 \right)^2}-...-\frac{1}{(mn)^2}$$ The second derivative is always negative.
That means the derivative is a decreasing function.
Also, $f'(m)\rightarrow \log{(1-p)}<0$ as $m\rightarrow \infty$ because there are a fixed, finite number of other terms in $f'(m)$ and each of those other terms converge to 0.
If the derivative at the lowest possible $m$, that is at $\max{X_i}$, is less than or equal to 0, then the mle is $\max{X_i}$. On the other hand, if $f'(\max{X_i})>0$, then let $M$ be the largest integer such that $f'(M)>0$. The mle is either $M$ or $M+1$. Just plug both of them into the likelihood function to see which of them makes $f$ bigger.

Note that
a) $\left(\max{X_i},\sum{X_i} \right)$ is a sufficient statistic
b) $$f'(m)=\log{(1-p)}+\sum_{j=1}^{\sum{X_i}}\frac{1}{mn-\sum{X_i}+j}\approx \log{(1-p)}+\int_{mn-\sum{X_i}+0.5}^{mn+0.5}x^{-1}dx$$ $$=\log{(1-p)}+\log{(mn+0.5)}-\log{(mn-\sum{X_i}+0.5)}$$ and the approximate solution to $f'(m)=0$ is therefore $m\approx \frac{\sum{X_i}-0.5 p}{n p}$.

Here are some examples from simulated data using the following R program:

mlem=function(X,p) {
  sX=sum(X)
  mtry=mX=max(X)
  if ((log(1-p)+sum(1/c((mtry*n+sX+1):(mtry*n))))<0) {
    mle=mtry
    mtry=mtry+1
    } else {
    while ((log(1-p)+sum(1/c((mtry*n-sX+1):(mtry*n))))>0) mtry=mtry+1
    if (dbinom(sX,n*(mtry-1),p)>dbinom(sX,n*mtry,p)) mle=mtry-1 else mle=mtry
  }
#return the mle, max Xi, Max integer M with f'(M)>0, approximate solution
  return(c(mle,mX,mtry-1,round((sX-0.5*p)/(length(X)*p))))
}


set.seed(123)
for (i in 1:3) {
  p=runif(1)
  m=30+round(20*runif(1))
  n=10+round(10*runif(1))
  X=rbinom(n,m,p)
  mle=mlem(X,p)
  print(c(m,n,p))
  print(X)
  print(mlem(X,p))
}

The first line shows the true values of the parameters used to simulate the data. The second line shows the simulated values of $X_i$:

Example 1.
m=46 n=14 p=0.2875775
17 18 8 13 17 14 13 19 13 15 14 9 17 11
mle= 49, M= 49

Example 2.
m=37 n=20 p=0.04205953
3 2 2 5 2 2 2 2 1 0 4 3 2 2 0 1 2 1 1 1
mle=45, M=45

Example 3.
m=38 n=14 p=0.1428
5 3 3 4 5 4 8 2 5 7 3 6 4 3
mle=31, M=30

These examples show that the mle is not always equal to $M$ or always equal to $M+1$ where $M$ is the largest integer with $f'$ positive; it can be either of those. In all three of those examples, the mle was the nearest integer to $\frac{\sum{X_i}-0.5 p}{n p}$.

Sometimes, the mle is equal to $\max{X_i}$ as in this example simulated from data with $M=10$, $N=20$, and $p=0.8$:
9 7 8 6 6 10 8 6 8 8 6 8 8 8 10 6 9 10 9 6
Here, $f'(10)\approx -0.104$ and so the mle is $10$. Again in this example, it turns out that the mle is the nearest integer to $\frac{\sum{X_i}-0.5 p}{n p}$.

Summary:
If $f'(\max{X_i})<0$, then the mle is $\max{X_i}$.
Otherwise, find the largest integer $M$ such that $f'(M)>0$ using the formula above. $f'(x)$ is a decreasing function that goes to $-\infty$ as $x\rightarrow \infty$, so it is guaranteed that such an $M$ can be found. That $M$ will be close to $\frac{\sum{X_i}-0.5 p}{n p}$. Either $M$ or $M+1$ will be the mle.

$\endgroup$
1
$\begingroup$

I don't see how to do this with only one observation, if $p$ is unknown. If $p$ is known, here are some clues. [More generally, this is a much-studied problem; perhaps see this paper and its references.]

If $p = 0.3$ and your observation is $X = 12,$ then the method of moments estimator is $X/p = 40.$

Perhaps it is reasonable to guess that the MLE will be approximately the same.

Here is a graphical solution, using R:

m = 1:150; p = .3; x = 12
like = dbinom(x, m, p)
plot(m, like, type="l", lwd=2)
 abline(h=0, col="green2")
mle = mean(m[like==max(like)]);  mle
[1] 39.5
 abline(v=mle, col="red")

enter image description here

$\endgroup$
1
  • $\begingroup$ in OP notation $p$ is $n$ (number of trials). $\endgroup$ Mar 1, 2021 at 20:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.