6
$\begingroup$

Let $X=\{x_1,x_2,...x_m\}$ be a random variable with density $f(x)=cx^{2n}$ if $-1\lt x\lt 1$, $0$ otherwise, and $n \in \{ 1,2,3,4\}$. Find the MLE of $n$.

I tried to find the Likelihood function : $(n+1/2)^m\prod_{i=1}^mx_i^{2n}$ but I'm not sure. Suggestions? Thanks

$\endgroup$
  • 1
    $\begingroup$ Given the likelihood function, all that's left to do is to maximize it with respect to n. Do you know how to do that? Hint: What are the possible values of n which have to be considered when maximizing the likelihood function? $\endgroup$ – Mark L. Stone May 7 '16 at 22:03
  • $\begingroup$ yes, Mark L. Stone, but is correct my likelihood function? $\endgroup$ – albert May 8 '16 at 15:43
4
$\begingroup$

This problem is simpler than it might look: although it might get confusing when one tries to apply routine Calculus methods, it is easy when worked from general principles.


By definition, the likelihood $\mathcal L$ is the probability of the data. Since the data are (implicitly) assumed independent, this is the product of the individual probability densities, each equal to $(n+1/2)(x_i^2)^n$. Consequently, as shown in the question,

$$\mathcal{L}(n) = \prod_{i=1}^m \left((n+1/2)x_i^{2n}\right) = (n+1/2)^m \left(\prod_{i=1}^m x_i^2\right)^n.$$

The only part depending on the data is that product on the right. Since it is almost surely positive for any of the $n$ under consideration, we may write it in terms of its logarithm. A convenient multiple of that is the mean log of the squared data, $$q=\frac{1}{m}\log\prod_{i=1}^m x_i^2.$$ Note that $q$ depends on the data, but not on the unknown parameter $n$: it will be the test statistic on which the Maximum Likelihood estimate is based.

Indeed, let's take logs of both sides, obtaining

$$\log\mathcal{L}(n) = m\log(n+1/2) + mnq.$$

Observe that since $|x_i| \lt 1$ for all $i$, $q \lt 0$. To maximize $\mathcal{L}$ we merely have to select the largest of the four values with $n=1,2,3,4$: namely,

$$\eqalign{\log \mathcal{L}(1) &= m\log(3/2) + mq,\\ \log \mathcal{L}(2) &= m\log(5/2) + 2mq,\\ \log \mathcal{L}(3) &= m\log(7/2) + 3mq,\\ \log \mathcal{L}(4) &= m\log(9/2) + 4mq.}$$

From row to row the $m\log(n+1/2)$ terms increase, but they do so more and more slowly; yet the $nmq$ terms decrease the value of $\log\mathcal L(n)$ at a constant rate of $m|q|$. Thus, the largest value of $(1/m)\log\mathcal L$--which corresponds to the largest value of $\mathcal L$ itself--is attained at the point where, in scanning these four values, $|q|$ first exceeds $$\log((n+1)+1/2) - \log(n+1/2) = \log\left(\frac{2n+3}{2n+1}\right).$$

This leads to an extremely simple procedure that can be specified even before collecting the data. Namely,

The increasing sequence $\log(5/3) \gt \log(7/5) \gt \log(9/7)$ partitions the real numbers into four intervals $$I_1=(-\infty, \log\frac{9}{7}],\ I_2=[\log\frac{9}{7}, \log\frac{7}{5}],\ I_3=[\log\frac{7}{5}, \log\frac{5}{3}],\ I_4=[\log\frac{5}{3},\infty).$$ If $|q|\in I_j$, pick $\hat n=j$ as the maximum likelihood estimator.

If $|q|$ lies at the common boundary of two intervals, they give two MLEs for $n$.


For each $n=1,2,3,4$, this procedure was applied to $1000$ independent samples of size $10$, then $1000$ more independent samples of size $500$. The tabulation of the estimates shows decent correspondence between $n$ and $\hat n$ for the small samples--which is as much as one might hope--and near perfect correspondence for the large samples--which is where the Maximum Likelihood method ought to perform well. Here are the results (which can be reproduced with the following R code):

$`m = 10`
    MLE=1 MLE=2 MLE=3 MLE=4
n=1   753   169    15     1
n=2   213   488   241    77
n=3    27   222   340   283
n=4     7   121   404   639

$`m = 500`
    MLE=1 MLE=2 MLE=3 MLE=4
n=1  1000     0     0     0
n=2     0   999     0     0
n=3     0     1   998     2
n=4     0     0     2   998

For instance, when the true parameter was $n=2$ and the sample size was only $m=10$, the top table shows the MLE was correct $488$ out of $1000$ times and was off by one (that is, estimating $n$ as either $1$ or $3$) $213+241$ other times. On the whole, the MLE appears a little biased towards the middle values for small $m$ and extremely accurate for large $m$.

#
# Generate random values.
#
rpow <- function(n, p) {
  q <- runif(n)
  i <- sample.int(2, n, replace=TRUE)
  ifelse(i==1, -1, +1) * q^(1/(2*p+1))
}
#
# Confirm `rpow` works as intended by matching histograms to the PDFs.
#
par(mfrow=c(1,4))
for (n in 1:4) {
  hist(rpow(1e4, n), main=paste("n =", n), freq=FALSE)
  curve(x^(2*n)*(n+1/2), add=TRUE, col="Red", lwd=2)
}
#
# Test the MLE.
#
MLE.pow <- function(x, breaks=log(c(5/3, 7/5, 9/7))) {
  # Given data `x`, return the MLE in {1,2,3,4}.
  m <- length(x)
  q <- sum(log(x^2)) / m
  sum(abs(q) <= breaks) + 1
}
m <- c(10,500) # Specify sample sizes to test
set.seed(17)   # Create a reproducible starting point
sim <- lapply(m, function(m) {
  x <- replicate(1e3, sapply(1:4, function(n) MLE.pow(rpow(m, n))))
  results <- apply(x, 1, tabulate, nbins=4) # Tabulate the MLEs
  rownames(results) <- paste0("n=", 1:4)
  colnames(results) <- paste0("MLE=", 1:4)
  results     
})
names(sim) <- paste("m =", m)
print(sim)
$\endgroup$
  • 1
    $\begingroup$ (+1) I really like it when a forgotten question suddenly becomes filled with useful content! $\endgroup$ – Alecos Papadopoulos Jun 15 '17 at 18:55
3
$\begingroup$

This has a few tricky points, so let's work it out.

First, we note that $2n$ is an even number, so the function that wants to be a density will be non-negative from that respect, as it should, even though the variable $X$ may take negative values.

Second we need to determine the value of $c$ so that we have a proper density. We require

$$\int_{-1}^1 cx^{2n}dx =1 \implies \frac{c}{2n+1}x^{2n+1} \Big|^1_{-1} =1$$

and since $2n+1$ is an odd number, we get

$$\frac{c}{2n+1}[1-(-1)] =1 \implies c= n+1/2$$

...as already noted in the comments.

Therefore

$$ f_X(x) = (n+1/2)x^{2n}, \;\;\; -1<x<1$$

and the likelihood function from an i.i.d. sample is ($I\{\}$ being the indicator function)

$$L = I_{\{-1<x<1\}}\cdot (n+1/2)^m\prod_{i=1}^mx_i^{2n} $$

We usually consider the log-likelihood for various beneficial reasons, but here, if we take logarithms we will be looking at $\ln x_i$ which is not defined if $x_i\leq 0$. But if one attempted to maximize the likelihood directly to avoid taking logs it would again hit on the same problem -the logarithm of the $x$-values would again appear if we considered the derivative with respect to $n$.

But since $x_i^{2n}$ is non-negative, we can write

$$L = I_{\{-1<x<1\}}\cdot (n+1/2)^m\prod_{i=1}^m|x_i|^{2n} $$

$$\implies\ln L = \ln \left(I_{\{-1<x<1\}}\right) + m\ln(n+1/2)+2n\sum_{i=1}^m\ln|x_i|$$

..which leaves us only with the problem of obtaining a realized value $x_i=0$ exactly. Well, from a theoretic point of view we invoke the zero-probability of a continuous r.v. taking a specific value, while from an applied point of view, if our sample contains an exact zero value, we can just discard it.

While this log-likelihood is strictly concave in $n$ and has a straightforward f.o.c

$$\frac{\partial \ln L}{\partial n} = \frac {m}{n+1/2} + 2\sum_{i=1}^m\ln|x_i|$$

$$\implies \hat n = \frac {m}{-2\sum_{i=1}^m\ln|x_i|} - \frac 12$$

...we have just ignored the restrictions on the values of $n$. This is a discrete-optimization problem. Nothing guarantees that the above expression will give an integer, or that it will fall in between $1$ and $4$.

After a little thinking, instead of considering derivatives etc as usual and dance around $\hat n$, it would be simpler to evaluate the log-likelihood at the four possible values of $n$ and see for which it attains the highest value.

$\endgroup$

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.