A function that produces a plot for each face of an empirical 2D
variogram based on residuals produced after the fitting of a model
using the function asreml.
It also adds envelopes to the plot by simulating data sets in parallel
from a multivariate normal distribution with expectation equal to the
fitted values obtained from the fixed and spline terms and variance
matrix equal to the fitted variance matrix
(Stefanova, Smith & Cullis, 2009). The plot is controlled by the
residual model, which must consist of two factors corresponding to
the two physical dimensions underlying the data. It can also have a third
term involving the at or dsum function that defines sections
of the data, such as experiments in different environments.
In this case, the two variogram faces are produced for each section.
# S3 method for asreml
variofaces(asreml.obj, means=NULL, V=NULL,
sections = NULL, row.factor = NULL, col.factor = NULL,
nsim=100, seed = NULL,
extra.matrix = NULL, ignore.terms = NULL, fixed.spline.terms = NULL,
bound.exclusions = c("F","B","S","C"), tolerance=1E-10,
units = "ignore", update = TRUE, trace = FALSE,
graphics.device=NULL, ncores = 2, ...)A list with the following components:
face1: a data.frame containing the variogram values on
which the plot for the first dimension is based.
face2: a data.frame containing the variogram values on
which the plot for the second dimension is based.
An asreml object from a call to asreml in which the
data argument has been set.
The vector of means to be used in generating simulated data sets. If
it is NULL, the fitted values based on object are used.
It must be the same length as the response variable for object.
The fitted variance matrix, i.e. having the appropriate pattern and values
given the model fitted to the observed data and the estimates of the
parameters obtained. If V is NULL then estimateV.asreml
is called to obtain it from asreml.obj
A single character string that specifies the name of the column
in the data.frame that contains the factor
that identifies different sections of the data to which separate spatial
models have been fitted.
A single character string nominating a factor
that indexes the rows of a grid that are one dimension of a
spatial correlation model. The factor must a column in
the data.frame stored in the asreml.obj.
A single character string nominating a factor
that indexes the columns of a grid that are one dimension of a
spatial correlation model. The factor must a column in
the data.frame stored in the asreml.obj.
The number of data sets to be simulated in obtaining the envelopes.
A single value, interpreted as an integer, that specifies the
starting value of the random number generator. The "L'Ecuyer-CMRG" random
generator is used and nextRNGStream is used to seed each core from the
original seed.
A matrix of order equal to the number of observations that is to
be added to the variance matrix, the latter based
on the information in asreml.obj. It is assumed that the sigma-parameterized
values of the variance parameter estimates, such as is given in the varcomp
component of summary.asreml, have been used in calculating
extra.matrix; the values in the vparameters component of
G.param and R.param may be either gamma- or sigma-parameterized.
The argument extra.matrix can be used in conjunction with
ignore.terms as a workaround to include components of the variance matrix
for variance functions that have not been implemented in estimateV.
A character giving terms from either the random or
residual models that are to be ignored in that their contributions to
the variance is not to be included in the estimated matrix. The term names are those
given in the vparameters component of the asreml object or the
varcomp component produced by summary.asreml, but only up to the
first exclamation mark (!). This can be used
in conjunction with estimateV.asreml as a workaround to include components
of the variance matrix for variance functions that have not been implemented
in estimateV.
A character vector giving one or more spline terms in the
random model that are regarded as fixed and so are to be ignored because
they are not regarded as contributing to the variance. The term names are those
given in the vparameters component of the asreml object or the
varcomp component produced by summary.asreml, but only up to the
first exclamation mark (!).
A character specifying one or more bound codes that
will result in a variance parameter in the random model being excluded
from contributing to the variance. If set to NULL then none will
be excluded.
The value such that eigenvalues less than it are considered to be zero.
A character indicating whether the BLUPs for units are
added to the residuals when this reserved factor is included in the
random model. Possible values are addtoresiduals and
ignore. If standardized conditional residuals are plotted and
the BLUPs for units are to be added then it is the standardized
BLUPs that are added.
If TRUE, and set.terms is NULL, then
newfit.asreml is called to fit the model to be tested,
using the values of the variance parameters stored in
the asreml.object, that is stored in asrtests.obj, as starting values.
If FALSE or set.terms is not NULL, then
newfit.asreml will not use the stored variance parameter
values as starting values when fitting the new model, the only
modifications being ((i) the model is fitted to simulated data and
(ii) those specified via ..., except that changes
cannot be made to any of the models.
If TRUE then partial iteration details are displayed when ASReml-R functions are invoked; if FALSE then no output is displayed.
A character specifying a graphics device for plotting.
The default is graphics.device = NULL, which will result
in plots being produced on the current graphics device. Setting it to
"windows", for example, will result in a windows graphics
device being opened.
A numeric specifying the number of cores to use in doing the
simulations. In choosing a value for ncores, it is necessary to
take into account other processes that are using parallel processing at
the same time.
Other arguments that are passed down to the function asreml. Changes
to the models are not allowed. Other changes are dangerous and generally
should be avoided.
Chris Brien
The residual model is scanned to ensure that it involves only two factors
not included in the at function, and to see if it has a third factor in
an at function. If so, the faces of the 2D variogram, each based on one
of the two non-at factors, are derived from the residuals in the
supplied asreml object using asreml.variogram, this yielding the observed
variogram faces. If aom was set to TRUE for the asreml
object, the standardized conditional residuals are used.
Then nsim data sets are generated by
adding the fitted.values, extracted from the asreml object,
to a vector of values randomly generated from a normal distribution with
expectation zero and variance matrix V. Each data set is analyzed
using the model in object and several sets are generated and analyzed
in parallel. The variogram values for the faces are
obtained using asreml.variogram stored. Note, if the analysis for a
data set does not converge in maxiter iterations, it is discarded and
a replacement data set generated. The value of maxiter can be specified
in the call to variofaces.asreml. Plots are produced for each face and
include the observed values and the 2.5%, 50% & 97.5% quantiles.
Stefanova, K. T., Smith, A. B. & Cullis, B. R. (2009) Enhanced diagnostics for the spatial analysis of field trials. Journal of Agricultural, Biological, and Environmental Statistics, 14, 392--410.
asremlPlus-package, asreml, newfit.asreml,
plotVariofaces.data.frame, simulate.asreml, set.seed.
if (FALSE) {
data(Wheat.dat)
current.asr <- asreml(yield ~ Rep + WithinColPairs + Variety,
random = ~ Row + Column + units,
residual = ~ ar1(Row):ar1(Column),
data=Wheat.dat)
current.asrt <- as.asrtests(current.asr, NULL, NULL)
current.asrt <- rmboundary.asrtests(current.asrt)
# Form variance matrix based on estimated variance parameters
s2 <- current.asr$sigma2
gamma.Row <- current.asr$gammas[1]
gamma.unit <- current.asr$gammas[2]
rho.r <- current.asr$gammas[4]
rho.c <- current.asr$gammas[5]
row.ar1 <- mat.ar1(order=10, rho=rho.r)
col.ar1 <- mat.ar1(order=15, rho=rho.c)
V <- gamma.Row * fac.sumop(Wheat.dat$Row) +
gamma.unit * diag(1, nrow=150, ncol=150) +
mat.dirprod(col.ar1, row.ar1)
V <- s2*V
#Produce variogram faces plot (Stefanaova et al, 2009)
variofaces(current.asr, V=V, ncores = parallel::detectCores())
}Run the code above in your browser using DataLab