In this answer, I would like to elaborate a little on Matthew's +1 answer regarding the GLS perspective on what the econometrics literature calls the random effects estimator.
GLS perspective
Consider the linear model
\begin{equation}
y_{it}=\alpha + X_{it}\beta+u_{it}\qquad i=1,\ldots,m,\quad t=1,\ldots,T
\end{equation}
If it held that $E(u_{it}\vert X_{it})=0$ we could simply estimate the model by pooled OLS, which amounts to ignoring the panel data structure and simply lump all $n=mT$ observations together.
We model the $u_{it}$ using the error-component model
\begin{equation}
u_{it}=\eta_i+\epsilon_{it}
\end{equation}
In matrix notation, the model can be written as
\begin{equation}
y=\alpha \iota_{mT}+X\beta+D\eta+\epsilon
\end{equation}
where $y$ and $\epsilon$ are $n$-vectors with typical elements $y_{it}$ and $\epsilon_{it}$, and $D$ is an $n \times m$ (one column per unit) matrix of dummy variables. $D$ is such that if a row corresponds to an observation belonging to unit $i$, then $D$ has a one in column $i$ and 0 else, $i=1,\ldots,m$.
We furthermore assume
$$
E(\epsilon\epsilon^\prime)=\sigma_\epsilon^2I
$$
The individual-specific effects $\eta$ must be independent of the $\epsilon_{it}$. The random-effects estimator, unlike the fixed effects (again, econometrics terminology) one, however additionally requires the stronger assumption that
\begin{equation}
E(\eta_i\vert X)=0
\end{equation}
Under this assumption, pooled OLS would be unbiased, but we can derive a GLS estimator. Assume that the $\eta_i$ are IID with mean zero and variance $\sigma^2_\eta$.
This assumption accounts for the term random effects. Assuming, moreover, that the two error components are independent, it is easy to see that
\begin{align*}
\operatorname{Var}(u_{it})&=\sigma^2_\eta+\sigma^2_\epsilon\\
\operatorname{Cov}(u_{it},u_{is})&=\sigma^2_\eta\\
\operatorname{Cov}(u_{it},u_{js})&=0\qquad\text{for all } i\neq j
\end{align*}
We then get the following $n\times n$ variance-covariance matrix $\Omega$:
$$
\Omega=
\begin{pmatrix}
\Sigma&O&\cdots&O\\
O&\Sigma&\cdots&O\\
\vdots&\vdots&&\vdots\\
O&O&\cdots&\Sigma
\end{pmatrix}
$$
Here,
$$
\Sigma=\sigma^2_\eta \iota\iota^\prime+\sigma^2_\epsilon I_T
$$
with $\iota$ a $T$-vector of ones. We may hence write
$$
\Omega=\sigma^2_\eta (I_m\otimes\iota\iota^\prime)+\sigma^2_\epsilon (I_m\otimes I_T)
$$
For the GLS estimator
$$
\hat\beta_{RE}=(X'\Omega^{-1}X)^{-1}X'\Omega^{-1}y
$$
we require $\Omega^{-1}$. To this end, let $J_T=\iota\iota^\prime$, $\bar J_T=J_T/T$ and $E_T=I_T-\bar J_T$. Then, write
$$
\Omega=T\sigma^2_\eta (I_m\otimes\bar J_T)+\sigma^2_\epsilon (I_m\otimes E_T)+\sigma^2_\epsilon (I_m\otimes \bar J_T)
$$
or, collecting terms with the same matrices,
$$
\Omega=(T\sigma^2_\eta+\sigma^2_\epsilon) (I_m\otimes\bar J_T)+\sigma^2_\epsilon (I_m\otimes E_T)
$$
Idempotency of $P=I_m\otimes\bar J_T$ and $Q=I_m\otimes E_T$ then allows us to show that
$$\Omega^{-1}=\frac{1}{\sigma^2_1}P+\frac{1}{\sigma^2_\epsilon}Q= -\frac{\sigma^2_\eta}{\sigma^2_1\sigma^2_\epsilon}(I_m\otimes\iota\iota^\prime) + \frac{1}{\sigma^2_\epsilon}(I_m\otimes I_T),$$
where $\sigma^2_1=T\sigma^2_\eta+\sigma^2_\epsilon$.
Gauss-Markov logic then explains why the random effects estimator may be useful, as it is a more efficient estimator than pooled OLS or fixed effects under the given assumptions (provided, which is a very big if in many panel data applications, that the $\eta_i$ are indeed uncorrelated with the regressors). In short, GLS is more efficient because the error covariance matrix is not homoskedastic in this model.
One can show that GLS estimate can be obtained by running OLS on the partially demeaned data:
$$(y_{it}-\theta \bar y_{i\cdot}) = (X_{it} - \theta \bar X_{i\cdot})\beta + (u_{it} - \theta u_{i\cdot}),$$
where $\theta = 1-\sigma_\eta/\sigma_1$. For $\theta=1$ one gets the fixed effect ("within") estimator. For $\theta\to -\infty$ one gets the "between" estimator. The GLS estimator is a weighted average between the two. (For $\theta=0$ one gets the pooled OLS estimator.)
Feasible GLS
To make an FGLS approach practical, we require estimators of $\sigma^2_1$ and $\sigma^2_\epsilon$. Baltagi, Econometric Analysis of Panel Data, p. 16 (quoting from the 3rd edition), discusses the following options on how to proceed.
Assume first we observe $u_{it}$. Then,
$$\hat\sigma^2_1=T\frac{1}{m}\sum_{i=1}^m\bar{u}_{i\cdot}^2$$
and
$$\hat\sigma^2_\epsilon=\frac{1}{m(T-1)}\sum_{i=1}^m\sum_{t=1}^T\left(u_{it}-\frac{1}{m}\sum_{i=1}^m\bar{u}_{i\cdot}\right)^2$$
would be good estimators of their parameters, with $\bar{u}_{i\cdot}$ the time-average corresponding to the obseravations of unit $i$.
The Wallace and Hussein (1969) approach consists of replacing $u$ with residuals of a pooled OLS regression (which, after all, still is unbiased and consistent under the present assumptions).
The Amemiya (1971) approach suggests using FE (or LSDV) residuals instead. As a computational matter, we impose the restriction that $\sum_i\eta_i=0$ to circumvent the dummy variable trap so as to be able to get $\hat\alpha=\bar y_{\cdot\cdot}-\bar X_{\cdot\cdot}'\hat\beta_{FE}$ with $\cdot\cdot$ denoting grand averages over $i$ and $t$ for the LSDV residuals $\hat u=y-\hat\alpha-X\hat\beta_{FE}$.
The default Swamy and Arora (1972) approach estimates
$$
\hat\sigma^2_\epsilon=[y'Q(I-X(X'QX)^{-1}X'Q)y]/[m(T-1)-K]
$$
and
$$
\hat\sigma^2_1=[y'P(I-Z(Z'PX)^{-1}Z'P)y]/[m-K-1]
$$
Here, $Z=(\iota_{mT}\quad X)$.
The Nerlove (1971) approach estimates $\sigma_\eta^2$ from $\sum_{i=1}^m(\hat\eta_i-\bar{\hat\eta})^2/(m-1)$ where the $\hat\eta_i$ are dummies from a fixed effects regression and $\hat\sigma^2_\epsilon$ is estimated from the within residual sums of squares from this regression, with $mT$ in the denominator.
I am also very surprised that these make such a big difference as shown by Randel's calculations!
EDIT:
Regarding the differences, the estimates of the error components may be retrived in the plm package, and indeed return vastly different results, explaining the difference in the point estimates for $\beta$ (as per @Randel's answer, amemiya throws an error that I did not attempt to fix):
> ercomp(stackY~stackX, data = paneldata, method = "walhus")
var std.dev share
idiosyncratic 21.0726 4.5905 0.981
individual 0.4071 0.6380 0.019
theta: 0.06933
> ercomp(stackY~stackX, data = paneldata, method = "swar")
var std.dev share
idiosyncratic 0.6437 0.8023 0.229
individual 2.1732 1.4742 0.771
theta: 0.811
> ercomp(stackY~stackX, data = paneldata, method = "nerlove")
var std.dev share
idiosyncratic 0.5565 0.7460 0.002
individual 342.2514 18.5000 0.998
theta: 0.9857
I suspect that the estimators of the error components are also not consistent in my example in the sister thread where I aim to demonstrate differences between FE and RE using data where the individual effects and $X$ are correlated. (In fact, they cannot be, because they ultimately drive away the RE estimate from the FE estimate as per the fact that RE is a weighted average of FE and between estimation with weights determined by the error component estimates. So, if RE is not consistent, that must ultimately be due to these estimates.)
If you replace the "offending" feature of that example,
alpha = runif(n,seq(0,step*n,by=step),seq(step,step*n+step,by=step))
by simply, say,
alpha = runif(n)
so random effects that are uncorrelated with $X$, you get RE point estimates for $\beta$ very close to the true value $\beta=-1$ for all variantes of estimating the error components.
References
Amemiya, T., 1971, The estimation of the variances in a variance-components model, International Economic Review 12, 1–13.
Baltagi, B. H., Econometric Analysis of Panel Data, Wiley.
Nerlove, M., 1971a, Further evidence on the estimation of dynamic economic relations from a time-series of cross-sections, Econometrica 39, 359–382.
Swamy, P.A.V.B. and S.S. Arora, 1972, The exact finite sample properties of the estimators of coefficients in the error components regression models, Econometrica 40, 261–275.
Wallace, T.D. and A. Hussain, 1969, The use of error components models in combining cross-section and time-series data, Econometrica 37, 55–72.
