Application of Maximum Likelihood estimation (MLE) to the step of Feasible Generalized Least Square (FGLS) I have the following regression
$$y = X\beta +u$$
where $y$ and $u$ are $(n\times 1)$ and $X$ is a fixed $(n \times k)$ matrix with full column rank and $\beta$ is an unknown $(k\times 1)$ vector of parameters.
$E(u)=0$ and $var(u)=\sigma^2 V$ with $(n\times n) $ symmetric and positive definite non-diagonal matrix $V$.
Let $\hat {\beta_{OLS}}= (X'X)^{-1} X'y$
Suppose $$u_t =\phi u_{t-1} +v_t $$ where $v_t$ is a white noise process with variance $\sigma_v^2$ and $\phi$ is an unknown parameter with $|\phi|<1$.
Describe how to get the feasible generalized least square (FGLS) estimator for $\beta$.
In the FGLS estimation procedure, we should use the maximized likelihood estimator $\hat {\phi_{MLE}}$, which should be obtained through the conditional log-likelihood  function set conditional on the first observation.

What I did
I obtain the variance-covariance matrix as follows:

Since $\phi$ is unknown, the matrix $V$ is also unknown. Therefore, we need to use FGLS as follows:
Step 1: Estimate the generalized regression $y=X\beta + u$ by OLS and obtain $\hat {u}$.
Step 2: Construct the auxiliary regression $\hat {u_t}= \phi \hat {u_{t-1}}+ v_t$ in order to obtain $\hat {\phi_{MLE}}$ and $\hat {v_{MLE}}$.
For that, I write down a likelihood function
$$v_t \sim N(0, \sigma_v^2)$$
$\hat {u_t} \sim N(0, \sigma_v^2)$
$I(\theta, u)=- \frac{n}{2} log(2\pi) - \frac{n}{2} log(\sigma_v^2) - \frac{n}{2\sigma_v^2} (\hat {u_t}- \phi \hat {u_{t-1}})'(\hat {u_t} - \phi \hat {u_{t-1}}) $
where $\theta = [\phi \sigma_v^2]$
FOC:
$$\frac{\partial I(\theta, u)}{\partial \phi} = 0$$
$$\hat {\phi_{MLE}} = (\hat {u_{t-1}}' \hat {u_{t-1}})\hat {u_{t-1}}' \hat {u_{t}}$$
Step 3:
By using $\hat {\phi_{MLE}}$ and $\hat {v_{MLE}}$, I will obtain the matrix $V$ and calculate $\hat {\beta_{FGLS}} = (X'V^{-1}X)^{-1}X'V^{-1}y$

I couldn’t do the step 2. The question states the conditional log-likelihood function, however I cannot write down this log-likelihood function as conditional function to solve for $\hat {\phi_{MLE}}$ and $\hat {v_{MLE}}$.
Please help me how to calculate the step 2. Thank you for your all kinds of helps and suggestions.

It would appear that Feasible Generalized Least Square (FGLS) should be used when the covariance matrix of the errors $cov(u) = \sigma^2V$ has a completely unknown form. However, in this problem, you have very specific prior knowledge of its form: it's the autocorrelation matrix of a 1st order autoregressive process, AR(1), which is a Toeplitz matrix. You even wrote this matrix explicitly.
Therefore, I feel the most appropriate and exact solution is to "directly" maximize the conditional log likelihood (conditioned on $\phi$), equivalently minimize the negative log likelihood conditioned on $\phi$, and then maximize over $\phi$:
For each specific value of $\phi$, the log likelihood of $u$ conditioned on $\phi$ is a multivariate Gaussian of $u$ ($u_{\phi} = y - X\beta_{\phi}$), but it's necessary to keep the $det(V_{\phi})$ term, as explained below. For a specific value of $\phi$, i.e., conditioned on $\phi$, this is maximized by the Generalized Least Squares (GLS): $\beta_{\phi} = (X'V^{-1}_{\phi}X)^{-1}X'V^{-1}_{\phi}y$.
But in order to maximize the overall unconditional (log) likelihood, it is necessary to compare likelihood values for different values of $\phi$, which requires to keep the $det(V_{\phi})$ term in the expression of the log likelihood.
Since the value of the unknown $\phi$ is in the bounded 1D range $-1 < \phi < 1$, it is a simple matter to maximize the unconditional likelihood over this range by a numerical 1D search (for example by a Fibonacci search for the maximum, etc.), or analytically by computing the $det(V_{\phi})$ as a function of $\phi$ and differentiating by $\phi$. (The analytical expression for $det(V_{\phi})$ as a function of $\phi$ can be obtained "easily", since $V$ has a very specific "easy" Toeplitz form.)
