Put and assume . Then
Thus generally has a periodic mean, and generally has a periodic variance. The edge cases allowed by matter: is stationary when , because the mean is then constant, and is stationary for or . In the latter case multiplying iid centered Gaussian white noise by leaves its iid distribution unchanged. For every , but , so is nonstationary. For , is nonstationary because is not constant.
Use the seasonal difference operator . Since ,
This is a stationary moving-average model of iid Gaussian white noise. Its covariance is at lag , at lags , and zero otherwise.
For the variance-modulated process the same operation gives
Its variance is and its covariance at lag is . Therefore it remains nonstationary for , and the operation doubles the marginal variance at nonzero seasons. Seasonal differencing does not remove periodic variance: it removes a periodic mean, but the periodic variance generally remains. Practical alternatives are a periodic model or seasonal variance standardization; at seasons with the observations are deterministically zero, so division by is not possible there. For the already stationary special cases , filtering preserves stationarity.
A periodically correlated process with period has its mean and two-time covariance unchanged when both times are shifted by . For , the mean is zero and
Hence is a periodically correlated process. Its covariance period can be smaller than : when is even, , so is already a covariance period.
For a centered Gaussian random field, Wick theorem takes the classical form of the Isserlis theorem: every even correlator is a sum over pairings of two-point correlators, and every odd correlator vanishes. For example, with ,
A -point correlator has pairings. Equivalently, all connected correlation functions beyond order two vanish. The covariance completely determines the centered Gaussian statistics; higher ordinary correlators can be nonzero, but contain no independent connected information. For a field with nonzero mean, apply these statements to the centered field and then restore its mean.
There is a genuine syntax error in the printed command: random=list(...) must be a separate argument, following a comma rather than a plus. The intended call is
model2 <- gamm(npres ~ s(age, bs="cr") + s(wbc, bs="cr"), random=list(doctor=~age))
With the default Gaussian family, this generalized additive mixed model is
Here across doctors, is an estimated covariance matrix, and the Gaussian errors have variance and are independent of the random effects. The formula ~age includes a random intercept as well as a random slope; it does not force these two effects to be uncorrelated.
The random intercept allows different prescribing baselines, and the random slope allows doctor-specific linear modifications to the population age curve. Conditional on these effects the observations are independent under the model, whereas two patients of the same doctor have additional covariance . This accounts for clustered data and avoids treating doctor variation as independent patient-level noise. The population smooth and doctor-specific predictions are different targets. With only four doctors, estimates of the random-effect covariance can be imprecise; the Gaussian support issue for counts also remains. The new model accounts for doctor-specific baselines and age slopes.
Let have columns , and let have one column per incubator, with entry in the column for observation 's incubator and zeros elsewhere. Integrating the Gaussian random effects gives
Consequently its marginal log-likelihood is
The likelihood-ratio test of against uses
For each covariance, the optimizing coefficients are the generalized least squares estimator . Both fits should use the same marginal likelihood convention, for example fitting the alternative with update(fly.model, REML=FALSE) and the null as a Gaussian linear model. The displayed alternative REML criterion alone cannot supply .
The variance component is constrained to be nonnegative, and zero lies on the boundary. Thus the interior-parameter hypothesis of Wilks theorem fails; the usual approximation is inappropriate. A limit sometimes applies to a single variance component under additional asymptotic conditions, but four incubators do not give a persuasive large-number-of-groups approximation.
Use a parametric bootstrap under the null: estimate there; simulate independent Gaussian errors with that fitted variance at exactly the observed times and grouping labels; refit both models by maximum likelihood for each simulated sample; and calculate . The upper-tail proportion, conventionally , estimates the p-value. Boundary fits with estimated must remain in the simulation distribution. A restricted-likelihood ratio and a correspondingly calibrated bootstrap are another option when the fixed-effect design is identical, but one must not mix ordinary and restricted likelihoods in a single statistic.
Add a location-specific random intercept to the random-slope linear mixed model:
Take the location effects, incubator slopes and individual errors to be mutually independent, with the original variances and . The random intercept models variation already present at hatching because of the origin of the eggs. Because larvae are subsequently assigned to incubators, location and incubator are crossed random effects, not automatically nested groups. The corresponding code is
fly.model.loc <- lmer(size ~ hours + (0+hours|incubator) + (1|location))
The resulting covariance between two observations contains if they share a location and if they share an incubator. If the three named locations themselves are the only populations of interest, use a fixed location factor instead:
fly.model.loc.fixed <- lmer(size ~ hours + location + (0+hours|incubator))
Choosing between these interpretations depends on the inferential target; with only three locations, a location variance component will be weakly estimated.
For a weakly stationary process, every observation has the same finite expectation . Linearity of expectation, without any assumption of independence, gives
The variance is a sum of all pairwise covariances:
For a real-valued weakly stationary process, the autocovariance is even. There are pairs at lag zero and pairs in each direction at positive lag . Therefore
Positive autocovariances inflate the uncertainty relative to the independent-observation formula; negative terms may reduce it. An unbiased estimator need not have variance tending to zero: a process consisting of the same random variable at every time is weakly stationary, but averaging cannot reduce its variance.
A zero-mean autoregressive process of order one satisfies
where is white noise with variance . A causal time series uses present and past innovations only. For the series
has mean-square convergence because , and it supplies the causal weakly stationary process. Its variance is . For , only the common innovations contribute to the autocovariance:
Evenness then gives
When , the signs alternate; when , the process is white noise. With nondegenerate innovations, does not give a causal finite-variance stationary solution. These covariance calculations need uncorrelated innovations, not Gaussianity; Gaussian innovations additionally yield a Gaussian process.
There is a literal domain error in the printed formula. A semivariogram indexed by a signed displacement must satisfy , since reversing the two observations leaves the squared increment unchanged. For positive or , the displayed positive-lag value is nonzero whereas the displayed negative-lag value is zero. Thus the claim on all of is false as written.
Interpret the argument as a nonnegative distance, or use the even extension of the Gaussian semivariogram. For and , the intended model is
It admits the covariogram
The Gaussian kernel is a positive-definite kernel, as can be seen from the Fourier transform of a Gaussian: in one dimension
and the nonnegative spectral density makes every finite covariance quadratic form nonnegative. In higher dimensions use the product Gaussian density. Adding independent Gaussian white noise of variance adds the diagonal nugget covariance. Thus a zero-mean stationary Gaussian random field exists with this covariance, and gives the corrected semivariogram.
For nonzero structured variance, the nugget effect is the right limit at zero, the sill of a semivariogram is , and the range of a semivariogram is infinite if defined as the distance at which the sill is exactly reached. A 95-percent practical range solves :
The scale is sometimes called the range parameter; it is neither the exact range nor the 95-percent practical range. If , there is no nontrivial structured range.
The fact that a semivariogram does not determine stationarity is a further distinction: a semivariogram specifies increment variation, not a unique process or absolute covariance. Adding an independent random constant adds a constant to without changing . More strongly, subtracting from a stationary field preserves its semivariogram but generally makes its variance depend on . The corrected statement is that this semivariogram admits a second-order stationary model, with the displayed covariance using the usual convention at large distance; it does not force every compatible process to be stationary.
Write , let , let , and put . Assume the covariance matrix is invertible. Because the mean is known to be zero, every linear predictor is unbiased. Its mean squared prediction error is
Completing the square yields
The positive-definite matrix property makes the last term nonnegative, proving the simple kriging formulas
The residual is uncorrelated with . Since the joint law is multivariate normal, it is independent of , so the same predictor and variance describe the conditional distribution. Consequently a 95-percent prediction interval is
There is no constraint that the weights sum to one: that is required for an unknown constant mean in ordinary kriging, not for the stated known-zero-mean problem. If is singular, use a Moore-Penrose inverse with the corresponding covariance compatibility condition. At an already observed location prediction of that very same field value has zero error; predicting an independent noisy replicate is a different target. Estimated covariance parameters make this a plug-in interval and generally add uncertainty that the formula does not include.
The first binned empirical semivariogram rises strongly through the displayed distances, approximately quadratically, with no visible sill. The spatial data plot shows higher depth at larger . A stationary residual covariance should not be asked to explain that deterministic drift.
For a model with stationary zero-mean residual field, each raw squared difference has expectation
With linear drift , the added term is . Its average within distance bins can grow with separation, producing the shape seen in the first plot. The empirical semivariogram uses half the average squared difference over the pairs in each bin, so removing only a common intercept does not remove this drift contribution.
The formula depth ~ y estimates and removes a linear drift before computing a residual empirical semivariogram. The second graph levels off near with a fitted zero nugget and a finite scale, making the Gaussian residual covariance a plausible working model. Estimate the drift and fit spatial dependence to residuals rather than the raw depth differences. Estimated residuals are not the true errors, so their semivariogram has estimation effects; sparse long-distance bins are also noisy. These graphs support the model but do not establish stationarity or isotropy.
The subpart refers to svgm and svgm2, while the code calls the objects smvg and smvg2. Also the printed code omits the line that fits gauss.model2 before plotting or using it. The supplied output must be interpreted as the result of that missing residual-variogram fit, rather than as an executable complete script.
The universal kriging model has a linear drift in and a zero-mean Gaussian residual field:
The Gau range parameter in this software is , since its semivariogram uses . The output gives , , and , hence . The exact Gaussian range is infinite; the 95-percent practical range is approximately in the coordinate distance units.
The BLUE=TRUE option asks for the best linear unbiased estimator of the drift, rather than the full field kriging prediction. Thus the two printed fitted trends give at and at . Therefore
There is no coefficient in the drift. The reported variances and are the variances of the corresponding estimated trend values, not full-field prediction variances. They do not determine the slope variance without its covariance with the intercept. The interpolated map can still vary with through the correlated residual prediction, even though the fitted deterministic drift depends only on .
Let have rows and let be the spatial covariance matrix with . Conditional on chosen covariance parameters, the best linear unbiased estimator of the drift coefficients is the generalized least squares estimator
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
Substitute to obtain a profile log-likelihood, maximize over and with a positive-definite matrix , and recover the drift coefficients from the optimizing covariance. Restricted maximum likelihood adds the design determinant term and uses in place of 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. Estimate drift using covariance-weighted regression, and spatial dependence from drift-adjusted variation, with the two stages or likelihood explicitly distinguished.