= Solution
Let $X$ have rows $(1,y_i)$ and let $\Sigma_\vartheta$ be the spatial <covariance matrix> with $\vartheta=(\tau^2,\sigma^2,a)$. Conditional on chosen <covariance> parameters, the <best linear unbiased estimator> of the drift coefficients is the <generalized least squares> estimator
$$
\boxed{\widehat\beta(\vartheta)=(X^T\Sigma_\vartheta^{-1}X)^{-1}X^T\Sigma_\vartheta^{-1}z.}
$$
The printed procedure first estimates the linear drift from `depth ~ y`, constructs a binned residual <empirical semivariogram>, fits a Gaussian <semivariogram> to it, and then uses the fitted <covariance> in <generalized least squares> trend estimation and <universal kriging>. The missing fit could be supplied as
``
gauss.model2 <- fit.variogram(smvg2, vgm(0.5, "Gau", 1, 0))
``
where those arguments are initial values, not the reported final estimates. Such <semivariogram> fitting typically minimizes a weighted sum of squared discrepancies between binned empirical values and model values; it is not automatically a joint Gaussian likelihood maximization. One can iterate: form residuals using the current drift, refit the <covariance> model, recompute the <generalized least squares> drift, and repeat until changes are small. The displayed commands themselves show a single sequence, not evidence that such iteration occurred.
A likelihood-based alternative estimates the two blocks coherently from
$$
\ell(\beta,\vartheta)=-\frac12\{n\log(2\pi)+\log\det\Sigma_\vartheta+(z-X\beta)^T\Sigma_\vartheta^{-1}(z-X\beta)\}.
$$
Substitute $\widehat\beta(\vartheta)$ to obtain a profile <log-likelihood>, maximize over $\tau^2,\sigma^2\ge0$ and $a>0$ with a <positive-definite matrix> $\Sigma_\vartheta$, and recover the drift coefficients from the optimizing <covariance>. <Restricted maximum likelihood> adds the design determinant term $\log\det(X^T\Sigma_\vartheta^{-1}X)$ and uses $n-2$ in place of $n$ in the normalizing term, up to a fixed-design constant. It can reduce <bias> from estimating the drift before the <covariance>. A zero fitted nugget is a legitimate boundary estimate, not proof that measurement error is absent. \b[Estimate drift using covariance-weighted regression, and spatial dependence from drift-adjusted variation, with the two stages or likelihood explicitly distinguished.]
Back to article page
