R: Functions for generating a multiresolution, compactly...
Functions for generating a multiresolution, compactly supported basis,
multiresolution covariance functions and simulating from these processes.
These functions support the LKrig function. Their main
function is to create and evaluate radial basis functions of varying
support on a nested set of regular grids. This series of grids forms
a multiresolution basis. The Gaussian process model is an expansion
in these basis functions where the basis coefficients follow a Markov
random field model for each resolution level. This family of
functions generate the basis using sparse matrices, evaluate the
covariance function of the process and also simulate realizations of
the process. LKrig.cov.plot is a useful function to get a quick
plot of the covariance function implied by a LatticeKrig
If TRUE the conversion of the sparse matrix to
spam format will have at least one element in each row. If there are
no elements explicitly given in obj then an element with value
zero is added. This technical detail is needed to accommodate the spam
format for sparse matrices.
If passed the covariance matrix will be multiplied by this
vector or matrix.
The point in the spatial domain that is used to evaluate the
covariance function. The evaluation is done on x and y transects
through the spatial domain intersecting at center and finding
the covariance with respect to this point. If NULL defaults to the
center of the spatial domain.
A list giving the memory requirements for
a sparse cholesky decomposition. See chol.spam for details and
A list with components that give the information
describing a multiresolution basis with a Markov random field used for
the covariance of the basis coefficients. This list is created in
LKrig or by LKrigSetup and returned in the output
object. (See section on returned Value below for this list's
If TRUE returns the marginal variance. Currently not
Number of points to evaluate the covariance function along
each transect of the spatial domain.
An object returned by LKrig or a sparse matrix in
row/column format passed to LKrig.spind2spam.
A sparse matrix of basis functions (rows index points for
evaluation, columns index basis functions).
If TRUE B is returned instead of the precision matrix
A sparse (spam format) precision matrix.
The wX matrix constructed in the LKrig function.
A two column matrix of 2-dimension locations to evaluate
basis functions or the first set of locations to evaluate the
covariance function or the locations for the simulated process. Rows
index the different locations: to be precise x1[i,1:2] are the
"x" and "y" coordinates for the i th location.
Second set of locations to evaluate covariance function.
Limits in x coordinate for evaluating the covariance
model. Default is the spatial domain.
Limits in y coordinate for evaluating the covariance
model.Default is the spatial domain.
If TRUE intermediate steps and other debugging
information are printed.
The basis functions are two-dimensional radial basis
functions based on the compactly supported stationary covariance function
(Wendland covariance) and centered on
regular grid points with the scaling tied to the grid spacing.
For a basis at the coarsest level, the grid centers are generated by
expanding the two sequences
into a regular grid of center points. The same spacing delta is
used in both directions. The unnormalized basis functions are
evaluated at locations x1 by finding the pairwise, radial
distances among centers and x1, scaling by
grid.info$delta * overlap and then evaluating with the function
name passed as BasisFunction. By default this is the 2-d
Wendland covariance of order 2. Perhaps the most important point
about the LKrig.basis is that it is designed to return a matrix
of all basis functions as a sequence of points. There is no need to
have a function that evaluates basis functions individually. In R
code for a set of locations x1 and a rectangular spatial domain
with ranges xmin, xmax, ymin ,ymax:
Note that there will be nrow(centers) basis functions generated
where the precise number depends on the range of the domain and the
choice of delta. The basis functions are indexed by the columns in
PHI and this is a convention throughout this package. There
will be nrow(x1) rows in PHI as each basis function will be
evaluated at each 2-d location.
The basis functions are then normalized by scaling the basis functions
at each location so that resulting marginal variance of the process is
1. This is done to coax the covariance model closer to a stationary
representation. It is easiest to express this normalization by pseudo
If Q is the precision matrix of the basis coefficients then
in R/LatticeKrig code:
where Omega is the unnormalized covariance matrix of the basis
Although correct, the code above is not an efficient algorithm to
compute the unnormalized process variance. First the normalization
can be done level by level rather than dealing with the entire
multiresolution process at once. Also it is important to work with the
precision matrix rather than the covariance. The function
LKrig.normalize.basis.slow takes advantage of the sparsity of the
precision matrix for the coefficients.
LKrig.normalize.basis.fast is even a more efficient version when
a.wght is constant for a given level and takes advantage of the
Kronecker structure on the precision matrix at each level.
The precision matrix for the basis coefficients at each resolution has
the form t(B)%*% B. These matrices for the individual levels
are assembled by LKrig.precision as the block diagonals of a
larger precision matrix for the entire vector of coefficients. Note
these matrices are created in a sparse format. The specific
entries in B, the object created by LKrig.MRF.precision, are a
first order Markov random field: without edge adjustments the
diagonal elements have the value a.wght and the first order
neighbors have the value -1.
Below we give more details on how the weights are determined.
Following the notation in Lindgren and Rue a.wght= 4 + k2
with k2 greater than or equal to 0. Some schematics for filling
in the B matrix are given below (values are weights for the SAR on the
lattice with a period indicating zero weights).
Previous versions of LatticeKrig considered an edge correction to
reflect other boundary conditions. We have found these corrections to
be numerically unstable, however, and so prefer at this time of
writing adding a buffer of lattice points and using the uncorrected
weights described above.
LKrig.basis: A matrix with number of rows equal to the rows of
x1 and columns equal to the number of basis functions
(LKinfo$m). Attached to the matrix is an info attribute that
contains the list passed as LKinfo. Usually this value is in
spam sparse matrix format.
LKrig.precision: For return.B ==FALSE a sparse, square
matrix with dimensions of the number of basis functions. For
return.B == TRUE the "B" SAR matrix is returned. This is useful
for checking this function.
LKrig.cov: If C=NA a cross covariance matrix
with dimensions nrow(x1) and nrow(x2) is used. If C is
passed the result of multiplying the cross covariance matrix times
C is used.
LKrigMarginalVariance: Gives the marginal variance of the
LatticeKrig process at each level and at the locations in x1. Returned value is a matrix with
nlevel columns indexing the levels and the number of rows equal to nrow(x1).
If varLevels is a row of this matrix then sum( varLevels* LKinfo$alpha) is the marginal variance
for the full process when the different levels are weighted by alpha. This is weighted sum is the marginal
variance returned by LKrig.cov and marginal=TRUE ( Also assuming that LKinfo$rho.object
is NULL, which it usaully is.) .
LKrig.sim: A matrix with dimensions of nrow(x1) by
M. Each column are vectors of simulated values at the locations
LKrig.cov.plot: Evaluates the covariance specified in the list
LKinfo with respect to the point center along a transects in
the x and y directions intersecting this point. Note the rectangular
extent of the spatial domain is part of the grid information in
LKinfo. Returns components u, d and cov. Each of
these are two column matrices with the first column being results in
the x direction and second column in the y direction. d are the
distances of the points from the center and u are the actual x
or y spatial coordinates. cov are the values of the covariance
function. If normalize is TRUE these will in fact be the correlation
functions. To plot the returned list use
LKrig.quadraticform: Returns a vector that is
diag(t(PHI)%*% solve(Q) %*% PHI)) closely related to the marginal
variance of the process.
LKrig.normalize.basis.slow: A vector of variances
corresponding to the unnormalized process at the locations.
LKrig.spind2spam: This converts a matrix in spind sparse format to spam sparse format. Although useful for checking, LatticeKrig now uses the
the spam function directly to do the conversion.
If obj is a list in spind format then the