How to show that GLS estimator is consistent in regression model? $$\mathbf Y=\mathbf X\beta + \mathbf U$$
$$\mathbf E[\mathbf U\mathbf U']=\Omega $$

$$\hat \beta - \beta=(\mathbf X'\Omega^{-1}\mathbf X)^{-1}\mathbf X'\Omega^{-1}\mathbf U$$

let $\mathbf X^{*}=\Omega^{-0.5}\mathbf X$, and $\mathbf U^{*}=\Omega^{-0.5}\mathbf U$. it is obviously that $\mathbf u^{*}_i$,which is the i-th element of $\mathbf U^{*}$, is iid follows the standard normal distribution, if we assume that $\mathbf U$ follows normal distribution. But $x_i^{*'}$, the i-th row of $\mathbf X^{*}$ is not necessarily independent, for it is the linear combination of all rows in $\mathbf X$ (which is assumed iid for they stand for observations from the same distribution).

Thus the followings are not necessarily true under independent Law of Large number and Central limit theorem:


$plim\frac{\mathbf X'\Omega^{-1}\mathbf U}{n}=plim\frac{\mathbf X^{*'}\mathbf U^{*}}{n}= plim\frac{1}{n}\Sigma x_i^{*}u_i^{*} =0$
$\frac{\mathbf X'\Omega^{-1}\mathbf U}{\sqrt n}=\frac{\mathbf X^{*'}\mathbf U^{*}}{\sqrt n}= \frac{1}{\sqrt n}\Sigma x_i^{*}u_i^{*} $, converges to normal distribution


I understand that when $\mathbf \Omega $ is a diagonal matrix then $x_i^{*}$ is independent and it is not a problem, what if when $\mathbf \Omega $ is not a diagonal matrix? What kind of dependent LLN and CLT can be applied here? should we make other assumptions?
If the regressors are strictly exogenous to the error term (which is the initial standard assumption made), meaning, if the error vector is mean-independent of the regressor matrix

$$E(\mathbf u \mid \mathbf X) = \mathbf 0$$

then consistency follows immediately. Set for convenience $\mathbf \Omega ^{-1/2} \equiv \mathbf C$. Then 

$$\text{plim}\left(n^{-1} (\mathbf C\mathbf X)'\mathbf u \right)= \text{plim}\left(n^{-1}\mathbf X'\mathbf C'\mathbf u\right)$$

This is a $k \times 1$ matrix. What would a typical element will be, say the first?  

The matrix $\mathbf C'\mathbf u$ is a $n\times 1$ vector,

$$\mathbf C'\mathbf u = \left(\sum_{i=1}^nc_{1i}u_i,..., \sum_{i=1}^nc_{ni}u_i\right)'$$

Mutliplied by the first row of $\mathbf X' = (x_{11}\;...\;x_{1n})$ gives

$$\left [\mathbf X'\mathbf C'\mathbf u\right]_1 =x_{11}\cdot \sum_{i=1}^nc_{1i}u_i\;+...\;+ x_{1n}\cdot\sum_{i=1}^nc_{ni}u_i$$

Now strict exogeneity comes into play: it is stronger and so implies that 

$$E(\mathbf u \mid \mathbf X) = \mathbf 0 \implies E(x_{1i}u_j)=0\;\;\forall i,j$$ 

So the above sum of sums, scaled by $n^{-1}$ will tend to zero asymptotically, since all products in each sum will tend to a zero expected value. This proves consistency.

Note that strict exogeneity is really needed, for the case where the matrix $\Omega$ is a function of the regressors. In such a case, we would be looking at products of the form

$$x_{1i}\cdot h(\mathbf x) \cdot u_j$$

which would still go to zero in expected value terms, because mean-independence implies orthogonality between the error term and any function of the regressors.

If $\Omega$ does not depend on the regressors, then we could weaken strict exogeneity to the RHS condition, which is  stronger than "contemporaneous uncorrelatedness", the latter needed for OLS. 

I  leave the case for asymptotic normality to anyone interested.