How to find maximum likelihood estimates of an integer parameter? H.W. Question: 

$x_1,x_2,\ldots,x_n$ are independent Gaussian variables with mean $\mu$ and variance $\sigma^2$. Define $y = \sum_{n=1}^{N} x_n$ where $N$ is unknown. We are interested in estimation of $N$ from $y$.

a. Given  $\hat N_1 = y/\mu$ determine its bias and variance.
b. Given $\hat N_2 = y^2/\sigma^2$ determine its bias and variance.
Ignoring the requirement for $N$ to be an integer
c. Is there an efficient estimator (look at both $\mu = 0$ and $\mu \ne 0$)?
d. Find the maximum likelihood estimate of $N$ from $y$.
e. Find CRLB of $N$ from $y$.
f. Does the mean squared error of the estimators $\hat N_1,\hat N_2$ attain CRLB when $N\to \infty$?

If anyone could direct me to the solution of the following problem it would be great.
Thanks,
Nadav
 A: The method whuber has used in his excellent answer is a common optimisation "trick" that involves extending the likelihood function to allow real values of $N$, and then using the concavity of the log-likelihood to show that the discrete maximising value is one of the discrete values on either side of a continuous optima.  This is one commonly used method in discrete MLE problems involving a concave log-likelihood function.  Its value lies in the fact that it is usually possible to get a simple closed-form expression for the continuous optima.
For completeness, in this answer I will show you an alternative method, which uses discrete calculus using the forward-difference operator.  The log-likelihood function for this problem is the discrete function:
$$\ell_y(N) = -\frac{1}{2} \Bigg[ \ln (2 \pi) + \ln (\sigma^2) + \ln (N) + \frac{(y-N\mu)^2}{N\sigma^2} \Bigg]
\quad \quad \quad
\text{for } N \in \mathbb{N}.$$
The first forward-difference of the log-likelihood is:
$$\begin{equation} \begin{aligned}
\Delta \ell_y(N) 
&= -\frac{1}{2} \Bigg[ \ln (N+1) - \ln (N) + \frac{(y-N\mu - \mu)^2}{(N+1)\sigma^2} - \frac{(y-N\mu)^2}{N\sigma^2} \Bigg] \\[6pt]
&= -\frac{1}{2} \Bigg[ \ln \Big( \frac{N+1}{N} \Big) + \frac{N(y-N\mu - \mu)^2 - (N+1)(y-N\mu)^2}{N(N+1)\sigma^2} \Bigg] \\[6pt]
&= -\frac{1}{2} \Bigg[ \ln \Big( \frac{N+1}{N} \Big) + \frac{[N(y-N\mu)^2 -2N(y-N\mu) \mu + N \mu^2] - [N(y-N\mu)^2 + (y-N\mu)^2]}{N(N+1)\sigma^2} \Bigg] \\[6pt]
&= -\frac{1}{2} \Bigg[ \ln \Big( \frac{N+1}{N} \Big) - \frac{(y + N \mu)(y-N\mu) - N \mu^2}{N(N+1)\sigma^2} \Bigg]. \\[6pt]
\end{aligned} \end{equation}$$
With a bit of algebra, the second forward-difference can be shown to be:
$$\begin{equation} \begin{aligned}
\Delta^2 \ell_y(N) 
&= -\frac{1}{2} \Bigg[ \ln \Big( \frac{N+2}{N} \Big) + \frac{2 N (N+1) \mu^2 + 2(y + N \mu)(y-N\mu)}{N(N+1)(N+2)\sigma^2} \Bigg] < 0. \\[6pt]
\end{aligned} \end{equation}$$
This shows that the log-likelihood function is concave, so its smallest maximising point $\hat{N}$ will be:
$$\begin{equation} \begin{aligned}
\hat{N} 
&= \min \{ N \in \mathbb{N} | \Delta \ell_y(N) \leqslant 0 \} \\[6pt]
&= \min \Big\{ N \in \mathbb{N} \Big| \ln \Big( \frac{N+1}{N} \Big) \geqslant \frac{(y + N \mu)(y-N\mu) - N \mu^2}{N(N+1)\sigma^2} \Big\}.
\end{aligned} \end{equation}$$
(The next value will also be a maximising point if and only if $\Delta \ell_y(\hat{N}) = 0$.)  The MLE (either the smallest, or the whole set) can be programmed as a function via a simple while loop, and this should be able to give you the solution pretty quickly.  I will leave the programming part as an exercise.
A: You started well by writing down an expression for the likelihood.  It is simpler to recognize that $Y,$ being the sum of $N$ independent Normal$(\mu,\sigma^2)$ variables, has a Normal distribution with mean $N\mu$ and variance $N\sigma^2,$ whence its likelihood is
$$\mathcal{L}(y,N) = \frac{1}{\sqrt{2\pi N\sigma^2}} \exp\left(-\frac{(y-N\mu)^2}{2N\sigma^2}\right).$$
Let's work with its negative logarithm $\Lambda = -\log \mathcal{L},$ whose minima correspond to maxima of the likelihood:
$$2\Lambda(N) = \log(2\pi) + \log(\sigma^2) + \log(N) + \frac{(y-N\mu)^2}{N\sigma^2}.$$
We need to find all whole numbers that minimize this expression.  Pretend for a moment that $N$ could be any positive real number.  As such, $2\Lambda$ is a continuously differentiable function of $N$ with derivative
$$\frac{d}{dN} 2\Lambda(N) = \frac{1}{N} - \frac{(y-N\mu)^2}{\sigma^2N^2} - \frac{2\mu(y-N\mu)}{N\sigma^2}.$$
Equate this to zero to look for critical points, clear the denominators, and do a little algebra to simplify the result, giving
$$\mu^2 N^2 + \sigma^2 N -y^2 = 0\tag{1}$$
with a unique positive solution (when $\mu\ne 0$)
$$\hat N = \frac{1}{2\mu^2}\left(-\sigma^2 + \sqrt{\sigma^4 + 4\mu^2 y^2}\right).$$
It's straightforward to check that as $N$ approaches $0$ or grows large, $2\Lambda(N)$ grows large, so we know there's no global minimum near $N\approx 0$ nor near $N\approx \infty.$  That leaves just the one critical point we found, which therefore must be the global minimum.  Moreover, $2\Lambda$ must decrease as $\hat N$ is approached from below or above.  Thus, 

The global minima of $\Lambda$ must be among the two integers on either side of $\hat N.$

This gives an effective procedure to find the Maximum Likelihood estimator: it's either the floor or the ceiling of $\hat N$ (or, occasionally, both of them!), so compute $\hat N$ and simply choose which of these integers makes $2\Lambda$ smallest.
Let's pause to check that this result makes sense.  In two situations there is an intuitive solution:


*

*When $\mu$ is much greater than $\sigma$, $Y$ is going to be close to $\mu,$ whence a decent estimate of $N$ would simply be $|Y/\mu|.$  In such cases we may approximate the MLE by neglecting $\sigma^2,$ giving (as expected) $$\hat N = \frac{1}{2\mu^2}\left(-\sigma^2 + \sqrt{\sigma^4 + 4\mu^2 y^2}\right) \approx \frac{1}{2\mu^2}\sqrt{4\mu^2 y^2} = \left|\frac{y}{\mu}\right|.$$

*When $\sigma$ is much greater than $\mu,$ $Y$ could be spread all over the place, but on average $Y^2$ should be close to $\sigma^2,$ whence an intuitive estimate of $N$ would simply be $y^2/\sigma^2.$  Indeed, neglecting $\mu$ in equation $(1)$ gives the expected solution $$\hat N \approx \frac{y^2}{\sigma^2}.$$
In both cases, the MLE accords with intuition, indicating we have probably worked it out correctly.  The interesting situations, then, occur when $\mu$ and $\sigma$ are of comparable sizes.  Intuition may be of little help here.

To explore this further, I simulated three situations where $\sigma/\mu$ is $1/3,$ $1,$ or $3.$  It doesn't matter what $\mu$ is (so long as it is nonzero), so I took $\mu=1.$  In each situation I generated a random $Y$ for the cases $N=2,4,8,16,$ doing this independently five thousand times.
These histograms summarize the MLEs of $N$.  The vertical lines mark the true values of $N$.

On average, the MLE appears to be about right.  When $\sigma$ is relatively small, the MLE tends to be accurate: that's what the narrow histograms in the top row indicate.  When $\sigma \approx |\mu|,$ the MLE is rather uncertain.  When $\sigma \gg |\mu|,$ the MLE can often be $\hat N=1$ and sometimes can be several times $N$ (especially when $N$ is small).  These observations accord with what was predicted in the preceding intuitive analysis.

The key to the simulation is to implement the MLE.  It requires solving $(1)$ as well as evaluating $\Lambda$ for given values of $Y,$ $\mu,$ and $\sigma.$  The only new idea reflected here is checking the integers on either side of $\hat N.$  The last two lines of the function f carry out this calculation, with the help of lambda to evaluate the log likelihood.
lambda <- Vectorize(function(y, N, mu, sigma) {
  (log(N) + (y-mu*N)^2 / (N * sigma^2))/2
}, "N") # The negative log likelihood (without additive constant terms)

f <- function(y, mu, sigma) {
  if (mu==0) {
    N.hat <- y^2 / sigma^2
  } else {
    N.hat <- (sqrt(sigma^4 + 4*mu^2*y^2) - sigma^2) / (2*mu^2)
  }
  N.hat <- c(floor(N.hat), ceiling(N.hat))
  q <- lambda(y, N.hat, mu, sigma)
  N.hat[which.min(q)]
} # The ML estimator

A: Comment: Here is a brief simulation in R for $\mu = 50, \sigma = 3,$
which should be accurate to 2 or three places, approximating the
mean and SD of $Y.$ You should be able to find $E(Y)$ and $Var(Y)$
by elementary analytic methods as indicated in my earlier Comment.
If we had $N = 100$ then $E(\hat N)$ seems unbiased for $N.$
N = 100;  mu = 50;  sg = 3
y = replicate( 10^6, sum(rnorm(N, mu, sg))/mu )
mean(y);  sd(y)
[1] 99.99997
[1] 0.6001208
N.est = round(y);  mean(N.est);  sd(N.est)
[1] 99.9998
[1] 0.6649131

