$Y_{n\times 1}\sim N(\mu {\bf 1}_n, \sigma^2V)$, where $V = (1-\rho)I_{n\times n} - \rho J_{n\times n}$, $I$ is identity matrix and $J$ is the matrix with entries all 1, ${\bf 1}_n$ is a vector of one's, and $\mu$ is a scalar . Basically $V$ is compound symmetric. I am trying find the maximum likelihood estimates (MLEs) of $\mu$, $\sigma^2$ and $\rho$. In my attempt, I got MLEs for $\mu$ and $\sigma^2$ as: $$\widehat{\mu} = \frac{{\bf 1}_n^{T}\widehat{V}^{-1}Y}{{\bf 1}_n^{T}\widehat{V}^{-1}{\bf 1}_n},$$ $$\widehat{\sigma}^2 = (Y - \widehat{\mu}{\bf 1}_n)^{T}\widehat{V}^{-1}(Y - \widehat{\mu}{\bf 1}_n)/n,$$ where $\widehat{V}$ is same as $V$ by replacing $\rho = \widehat{\rho}$, and $\widehat{\rho}$ is the MLE of $\rho$.
First, is it possible to get an analytical solution for $\widehat{\rho}$? Afer replacing $\widehat{\mu}$ and $\widehat{\sigma}^2$ in the log likelihood, I tried to optimize with respect to $\rho$. But the result seems incorrect. Please suggest any method to find MLE of $\rho$.