Finite-Sample Inference After Feasible GLS

Why Finite-Sample Inference Needs Care

Feasible generalized least squares (FGLS) replaces an unknown error covariance matrix with an estimate. This makes GLS usable when the covariance structure is not known, but it introduces an additional source of uncertainty: the estimated weights are random. In a large sample, plug-in inference can often be justified by asymptotic theory. In a small sample, treating the estimated covariance matrix as if it were known can give a reference distribution that is too optimistic or otherwise poorly calibrated.

The central issue is not that a plug-in calculation is necessarily wrong. It is that its familiar finite-sample justification may no longer apply. With known covariance and normally distributed errors, GLS estimates can have exact normal or $t$ distributions under specified conditions. Once covariance parameters are estimated from the same data, that exact result generally does not carry over. The size of the discrepancy depends on the model, the sample size, and how much information the data provide about the covariance parameters.

Definition: In FGLS, the covariance matrix $\Omega(\gamma)$ depends on an unknown parameter vector $\gamma$, which is estimated by $\hat{\gamma}$. The estimator $$ \hat{\beta}_{FGLS} = \left(X'\Omega(\hat{\gamma})^{-1}X\right)^{-1} X'\Omega(\hat{\gamma})^{-1}y $$ uses the estimated covariance structure in place of the unknown one. A conventional plug-in covariance estimate substitutes $\hat{\gamma}$ into the GLS covariance formula and then proceeds as though that substitution had not added uncertainty.

A plug-in standard error is a useful estimate of variability, but it does not by itself establish an exact $t$ reference distribution. The asymptotic normality of FGLS and plug-in inference are appropriate large-sample foundations; finite-sample procedures need additional structure or approximation.

A Small-Sample Approximation for a Weighted Estimate

One transparent special case is a fixed weighted average of independent group means. It is useful because the contribution of each estimated variance can be seen separately. Suppose group $g$ has $n_g$ independent observations from a normal population with mean $\mu_g$ and variance $\sigma_g^2$. Assume the groups are independent, and let $c_g$ be fixed coefficients. The target and its estimator are $$ \theta=\sum_{g=1}^{G}c_g\mu_g, \qquad \hat{\theta}=\sum_{g=1}^{G}c_g\bar{y}_g. $$ The variance of $\hat{\theta}$ is $Q=\sum_g c_g^2\sigma_g^2/n_g$. Replacing each $\sigma_g^2$ by its group sample variance $s_g^2$ gives $\hat{Q}=\sum_g q_g$, where $q_g=c_g^2s_g^2/n_g$.

Under the stated normal sampling assumptions, the sample mean and sample variance in each group are independent, and $(n_g-1)s_g^2/\sigma_g^2$ has a chi-squared distribution with $n_g-1$ degrees of freedom. These facts motivate a Satterthwaite approximation: match the first two moments of the estimated variance to a scaled chi-squared variable, then use the resulting degrees of freedom for a $t$ reference distribution.

Derivation: Let $d_g=n_g-1$. Let $Q_g=c_g^2\sigma_g^2/n_g$. Since $q_g=c_g^2s_g^2/n_g$, $E(q_g)=Q_g$ and $\operatorname{Var}(q_g)=2Q_g^2/d_g$, and independent groups give $$ E(\hat{Q})=Q, \qquad \operatorname{Var}(\hat{Q})=2\sum_g\frac{Q_g^2}{d_g}. $$ Approximate $\hat{Q}$ by $Q\chi^2_{\nu}/\nu$, whose mean is $Q$ and variance is $2Q^2/\nu$. Equating this variance to its plug-in estimate $2\sum_g q_g^2/d_g$, and substituting $\hat{Q}$ for $Q$, gives $$ \hat{\nu}_{S} = \frac{\left(\sum_g q_g\right)^2} {\sum_g q_g^2/d_g}. $$ This moment matching does not make the approximation exact. It supplies an estimated degrees of freedom for a small-sample $t$ approximation when the independent normal-group setup is credible.

Worked Example: A Weighted Mean With Unequal Coefficients

Consider two independent groups with four observations each. All values below are hypothetical. In group 0, the observations are 9, 11, 10, and 10; in group 1, they are 12, 14, 13, and 13. Estimate $\theta=0.9\mu_0+0.1\mu_1$, and test $H_0:\theta=10$.

The group means are $\bar{y}_0=10$ and $\bar{y}_1=13$, so $\hat{\theta}=0.9(10)+0.1(13)=10.3$. In each group the squared deviations from the mean sum to 2. Thus $s_0^2=s_1^2=2/3$, using three degrees of freedom per group.

The estimated variance contributions are $$ q_0=\frac{0.9^2(2/3)}{4}=0.135, \qquad q_1=\frac{0.1^2(2/3)}{4}=0.001667, \qquad \hat{Q}=q_0+q_1=0.1367. $$ Therefore $SE=\sqrt{0.1366667}=0.3697$, and the observed statistic is $$ t=\frac{10.3-10}{0.369685}=0.8115. $$

The Satterthwaite degrees of freedom are $$ \hat{\nu}_{S} = \frac{(0.135+0.0016667)^2}{0.135^2/3+0.0016667^2/3} \approx 3.074. $$ Using a $t$ reference distribution with this approximate degrees of freedom gives a two-sided p-value of about 0.4752. The corresponding approximate 95% interval uses $t_{0.975,3.074}\approx3.14$: $$ 10.3\ \pm\ 3.14(0.3697) = [9.139,\ 11.461]. $$ Under the model and approximation, the interval is for the weighted mean in the original units of the observations. It includes the null value 10, consistent with the two-sided test at the 5% level.

For comparison, using the large-sample normal reference gives $0.8115$ a two-sided p-value of about 0.4171 and a 95% interval $10.3\pm1.96(0.3697)=[9.575,11.025]$. The difference illustrates how a small effective degrees of freedom can make inference more conservative. Neither calculation establishes that a general FGLS statistic has an exact $t$ distribution.

When the Weights Are Estimated

The fixed-coefficient calculation above isolates variance estimation while keeping the coefficients $c_g$ fixed. In FGLS, estimated covariance parameters can also change the weights used to calculate the coefficient estimate itself. The resulting estimator is a function of the observed outcomes through both the regression fit and the covariance estimate. Its finite-sample distribution need not be captured by simply substituting the estimated covariance matrix into the known-covariance formula.

A simple model makes the issue concrete. Suppose observations in two independent groups share a common mean $\beta$, but the group error variances can differ. If the variances were known, group $g$'s precision for estimating the common mean would be $n_g/\sigma_g^2$. FGLS instead uses $n_g/s_g^2$. When a group’s estimated variance is unusually small, that group receives more weight; with only a few observations, the weight itself can be highly variable.

Worked Example: Noisy Precision Weights in a Common-Mean Regression

Suppose the two groups have four observations each, and all observations are hypothetical. Group 0 has values 9, 10, 10, and 11. Group 1 has values 12, 12, 13, and 13. Assume independent normal errors within and across groups, a common mean $\beta$, and a separate variance for each group.

The sample means are $\bar{y}_0=10$ and $\bar{y}_1=12.5$. In group 0 the squared deviations sum to 2, giving $s_0^2=2/3$. In group 1 they sum to 1, giving $s_1^2=1/3$. The estimated group precisions are $4/(2/3)=6$ and $4/(1/3)=12$.

The FGLS estimate of the common mean is the precision-weighted average: $$ \hat{\beta}_{FGLS} = \frac{6(10)+12(12.5)}{6+12} = \frac{210}{18} = 11.667. $$ The conventional plug-in standard error is $$ SE_{\text{plug-in}} = \sqrt{\frac{1}{6+12}} = \sqrt{\frac{1}{18}} = 0.2357. $$

This standard error reflects the estimated group precisions as if they were fixed. But those precisions came from only four observations per group. For normally distributed observations, a sample mean and its sample variance are independent within a group, yet the variance estimate still fluctuates across repeated samples; the fitted estimate and its reported plug-in standard error therefore need not produce an exact finite-sample $t$ statistic. The displayed coefficient and standard error are a correct plug-in calculation, not an exact finite-sample confidence statement. A small-sample method should reflect the specific covariance-estimation procedure rather than assume that a generic degrees-of-freedom adjustment automatically applies.

Small-Sample Alternatives and Bootstrap Inference

There is no single exact small-sample correction for every FGLS model. The Satterthwaite approximation is particularly useful when uncertainty can be expressed as a sum of independent estimated variance components, as in the weighted-mean example. For a more complicated covariance model, a parametric bootstrap can reproduce the full estimation process: generate data under the null model, estimate the covariance parameters again in every simulated data set, refit FGLS, and recompute the test statistic.

Re-estimating the covariance parameters in each replicate is essential. Holding the original estimated weights fixed would mimic inference conditional on those weights, not the repeated-sampling procedure that produced them. A bootstrap still depends on the assumed data-generating model and is not automatically exact; it approximates the sampling distribution under that model.

Fit under the null.
Specify $H_0$ and estimate the nuisance parameters required by the null model, including covariance parameters.

Generate and refit.
Simulate a data set under the fitted null model. In that replicate, re-estimate the covariance parameters and compute the FGLS statistic using the same procedure as for the observed data.

Compare with the observed statistic.
For a two-sided test, count replicates with an absolute statistic at least as large as the observed absolute statistic. The proportion estimates the null tail probability.

Worked Example: A Null Bootstrap p-Value

Suppose a hypothetical FGLS analysis produces a two-sided observed statistic with absolute value $2.10$. In a parametric bootstrap, 999 data sets are generated under the null hypothesis, and the covariance model is re-estimated and the regression refitted for every data set. Suppose 63 simulated statistics have absolute value at least 2.10.

Using the plus-one estimate of the bootstrap p-value gives $$ \hat{p}_{\text{boot}} = \frac{63+1}{999+1} = 0.064. $$ The added one prevents a reported value of zero when no simulated statistic exceeds the observed one. For this illustrative run, the estimated two-sided p-value is 0.064, so it is not below a 5% significance level. This is a simulation-based approximation under the specified null model, not the probability that $H_0$ is true.

Common Pitfalls and Practical Interpretation

The most common mistake is to treat a plug-in covariance estimate as known merely because software reports a standard error and a test statistic. That reasoning confuses the estimated variance formula with the distribution of the statistic. The correction is to state what justifies the reference distribution: an exact result under a specified model, a stated small-sample approximation, or asymptotic normality.

Do not call a plug-in $t$ statistic exact without justification. Estimating covariance parameters generally changes finite-sample behavior.
Do not apply Satterthwaite degrees of freedom mechanically. The formula relies on a variance that is a sum of suitably independent estimated components; it is not a universal FGLS correction.
Do not bootstrap with the original weights frozen. Re-estimating covariance parameters in each replicate is part of reproducing FGLS uncertainty.
Do not interpret a p-value as a probability that the null is true. It is a tail probability for the statistic under the null model and the specified procedure.

In practice, report the covariance model, the method used to estimate its parameters, and the reference distribution or resampling procedure used for inference. When sample size is limited, comparing a conventional plug-in result with a defensible small-sample approximation or a carefully specified bootstrap can reveal whether the conclusion is sensitive to finite-sample uncertainty.

Check Your Understanding

Answer each question using the assumptions and calculations developed above.

Why does substituting an estimated covariance matrix into a GLS standard error not, by itself, establish an exact finite-sample $t$ distribution?
In the two-group weighted-mean example, what features of the data justify using the Satterthwaite approximation, and why is it still an approximation?
If the estimated variance contribution $q_0$ becomes much larger than the other contributions, what happens qualitatively to the effective degrees of freedom?
In a parametric bootstrap for FGLS, why must the covariance parameters be re-estimated in every simulated data set?
A two-sided bootstrap p-value is 0.064. State its interpretation and what it implies at a 5% significance level.