Package 'asremlPlus'

Description

Assists in automating the selection of terms to include in mixed models when 'asreml' is used to fit the models. Procedures are available for choosing models that conform to the hierarchy or marginality principle, for fitting and choosing between two-dimensional spatial models using correlation, natural cubic smoothing spline and P-spline models. A history of the fitting of a sequence of models is kept in a data frame. Also used to compute functions and contrasts of, to investigate differences between and to plot predictions obtained using any model fitting function. The content falls into the following natural groupings: (i) Data, (ii) Model modification functions, (iii) Model selection and description functions, (iv) Model diagnostics and simulation functions, (v) Prediction production and presentation functions, (vi) Response transformation functions, (vii) Object manipulation functions, and (viii) Miscellaneous functions (for further details see 'asremlPlus-package' in help). The 'asreml' package provides a computationally efficient algorithm for fitting a wide range of linear mixed models using Residual Maximum Likelihood. It is a commercial package and a license for it can be purchased from 'VSNi' <https://vsni.co.uk/> as 'asreml-R', who will supply a zip file for local installation/updating (see <https://asreml.kb.vsni.co.uk/>). It is not needed for functions that are methods for 'alldiffs' and 'data.frame' objects. The package 'asremPlus' can also be installed from <http://chris.brien.name/rpackages/>.

Version: 4.4.40

Date: 2024-10-27

Index

Note that many of the function below are S3 methods so that the suffix can be omitted. Of course, whether or not the suffix is omitted, the object supplied to the first argument must be of the class specified by the suffix. For example getFormulae.asreml is a getFormulae method for an asreml.object and so .asreml can be omitted and the object supplied to the first argument must be of class asreml.

The functions whose names end in 'alldiffs" utilize an alldiffs.object that stores: (i) a predictions.frame, being a data frame containing predicted values, variables indexing them and their standard errors and estimability status; the lower and upper limits of error intervals will be included when these are requested, (ii) optionally, square matrices containing all pairwise differences, the standard errors and p-values of the differences, and a data.frame containing LSD values and their summary statistics, (iii) optionally, the variance matrix of the predictions, and (iv) if the response was transformed for analysis, a data frame with backtransforms of the predicted values.

The functions whose names end in 'asrtests', which are most of the model functions, utilize an asrtests.object that stores: (i) the currently fitted model in asreml.obj, (ii) the table of test statistics for the fixed effects in wald.tab, and (iii) a data frame that contains a history of the changes made to the model in test.summary.

Author(s)

Chris Brien [aut, cre] (<https://orcid.org/0000-0003-0581-1817>)

Maintainer: Chris Brien <[email protected]>

References

Butler, D. G., Cullis, B. R., Gilmour, A. R., Gogel, B. J. and Thompson, R. (2023). ASReml-R Reference Manual Version 4.2. VSN International Ltd, https://asreml.kb.vsni.co.uk/

See Also

asreml

Examples

## Not run: 
## Analyse wheat dat using asreml and asremlPlus (see the WheatSpatial Vignette for details)
## Set up for analysis
library(dae)
library(asreml)
library(asremlPlus)
## use ?Wheat.dat for data set details
data(Wheat.dat)

# Add row and column covariates for the spatial modelling
tmp.dat <- within(Wheat.dat, 
                  {
                    cColumn <- dae::as.numfac(Column)
                    cColumn <- cColumn  - mean(unique(cColumn))
                    cRow <- dae::as.numfac(Row)
                    cRow <- cRow - mean(unique(cRow))
                  })

# Fit an initial model - Row and column random
current.asr <- do.call(asreml, 
                       list(yield ~ Rep + WithinColPairs + Variety, 
                            random = ~ Row + Column,
                            residual = ~ Row:Column,
                            data = tmp.dat))

# Intialize a model sequence by loading the current fit into an asrtests object
current.asrt <- as.asrtests(current.asr, NULL, NULL, IClikelihood = "full", 
                            label = "Initial model")

# Check for and remove any boundary terms and print a summary of the fit in the asrtests object 
current.asrt <- rmboundary(current.asrt)
print(current.asrt)


## Compare a series of information criteria to select a linear mixed model for the data

# Check the need for the term for within Column pairs (a post hoc factor)
current.asrt <- changeModelOnIC(current.asrt, dropFixed = "WithinColPairs", 
                                label = "Try dropping withinColPairs", IClikelihood = "full")
print(current.asrt)

# Fit an ar1 model for local spatial variation
spatial.ar1.asrt <- addSpatialModelOnIC(current.asrt, spatial.model = "corr", 
                                        row.covar = "cRow", col.covar = "cColumn", 
                                        row.factor = "Row", col.factor = "Column", 
                                        IClikelihood = "full")
spatial.ar1.asrt <- rmboundary(spatial.ar1.asrt)
infoCriteria(list(nonspatial = current.asrt$asreml.obj, 
                  ar1 = spatial.ar1.asrt$asreml.obj))
print(spatial.ar1.asrt)

# Choose a model for local spatial variation from several potential models
suppressWarnings(
  spatial.asrts <- chooseSpatialModelOnIC(current.asrt, 
                                          row.covar = "cRow", col.covar = "cColumn",
                                          row.factor = "Row", col.factor = "Column",
                                          dropRandom = "Row + Column",
                                          rotateX = TRUE, ngridangles = NULL, 
                                          asreml.option = "grp", return.asrts = "all"))
# Output the results
print(spatial.asrts$spatial.IC)
print(R2adj(spatial.asrts$asrts$TPNCSS$asreml.obj, include.which.random = ~ .))
print(spatial.asrts$best.spatial.mod)
print(spatial.asrts$asrts$TPNCSS)
printFormulae(spatial.asrts$asrts$TPNCSS$asreml.obj)

## Diagnosting checking using residual plots and variofaces

# Get current fitted asreml object and update to include standardized residuals
current.asr <- spatial.asrts$asrts$TPNCSS$asreml.obj
current.asr <- update(current.asr, aom=TRUE)
Wheat.dat$res <- residuals(current.asr, type = "stdCond")
Wheat.dat$fit <- fitted(current.asr)

# Do residuals-versus-fitted values plot
with(Wheat.dat, plot(fit, res))

# Plot variofaces
variofaces(current.asr, V=NULL, units="addtores", 
           maxiter=50, update = FALSE, 
           ncores = parallel::detectCores())

# Plot normal quantile plot
ggplot(data = Wheat.dat, mapping = aes(sample = res)) +
  stat_qq_band(bandType = "ts") + stat_qq_line() + stat_qq_point() +
  labs(x = "Theoretical Quantiles", y = "Sample Quantiles",
       title = "Normal probability plot") +
  theme(plot.title = element_text(size = 12, face = "bold")) + theme_bw()

## Prediction production and presentation

# Get Variety predictions and all pairwise prediction differences and p-values
Var.diffs <- predictPlus(classify = "Variety", 
                         asreml.obj=current.asr, 
                         error.intervals="halfLeast",
                         wald.tab=current.asrt$wald.tab, 
                         sortFactor = "Variety",
                         tables = "predictions")

# Plot the Variety predictions, with halfLSD intervals, and the p-values
plotPredictions(Var.diffs$predictions, 
                classify = "Variety", y = "predicted.value", 
                error.intervals = "half")
plotPvalues(Var.diffs)

## End(Not run)

Adds or recalculates the backtransforms component of an alldiffs.object.

Description

Given an alldiffs.object, adds or recalculate its backtransforms component. The values of transform.power, offset, scale and transform.function from the backtransforms component will be used, unless this component is NULL when the values supplied in the call will be used.

Usage

## S3 method for class 'alldiffs'
addBacktransforms(alldiffs.obj, 
                  transform.power = 1, offset = 0, scale = 1, 
                  transform.function =  "identity", ...)

Arguments

Value

An alldiffs.object with components predictions, vcov, differences, p.differences, sed, LSD and backtransforms.

The backtransforms component will have the attributes (i) LSDtype, LSDby and LSDstatistic added from the predictions component and (ii) transform.power, offset, scale, and link.

Author(s)

Chris Brien

See Also

asremlPlus-package, as.alldiffs, sort.alldiffs, subset.alldiffs, print.alldiffs,
renewClassify.alldiffs, redoErrorIntervals.alldiffs, plotPredictions.data.frame,
predictPlus.asreml, predictPresent.asreml

Examples

##Subset WaterRunoff data to reduce time to execute
data(WaterRunoff.dat)
tmp <- subset(WaterRunoff.dat, Date == "05-18" & Benches != "3")

##Use asreml to get predictions and associated statistics

## Not run: 
asreml.options(keep.order = TRUE) #required for asreml-R4 only
current.asr <- asreml(fixed = log.Turbidity ~ Benches + (Sources * (Type + Species)), 
                      random = ~ Benches:MainPlots,
                      keep.order=TRUE, data= tmp)
current.asrt <- as.asrtests(current.asr, NULL, NULL)
TS.diffs <- predictPlus(classify = "Sources:Type", 
                        asreml.obj = current.asr, 
                        wald.tab = current.asrt$wald.tab, 
                        present = c("Sources", "Type", "Species"))

## End(Not run)

##Use lmeTest and emmmeans to get predictions and associated statistics

if (requireNamespace("lmerTest", quietly = TRUE) && 
    requireNamespace("emmeans", quietly = TRUE))
{
  m1.lmer <- lmerTest::lmer(log.Turbidity ~ Benches + (Sources * (Type + Species)) + 
                              (1|Benches:MainPlots),
                            data=tmp)
  TS.emm <- emmeans::emmeans(m1.lmer, specs = ~ Sources:Species)
  TS.preds <- summary(TS.emm)
  den.df <- min(TS.preds$df, na.rm = TRUE)
  ## Modify TS.preds to be compatible with a predictions.frame
  TS.preds <- as.predictions.frame(TS.preds, predictions = "emmean", 
                                   se = "SE", interval.type = "CI", 
                                   interval.names = c("lower.CL", "upper.CL"))
  
  ## Form an all.diffs object and check its validity
  TS.vcov <- vcov(TS.emm)
  TS.diffs <- allDifferences(predictions = TS.preds, classify = "Sources:Species", 
                             vcov = TS.vcov, tdf = den.df)
  validAlldiffs(TS.diffs)
}  

## Recalculate the back-transforms of the predictions obtained using asreml or lmerTest
if (exists("TS.diffs"))
{
  TS.diffs <- addBacktransforms.alldiffs(TS.diffs, transform.power = 0)
}

Description

Adds either a correlation, two-dimensional tensor-product natural cubic smoothing spline (TPNCSS), or a two-dimensional tensor-product penalized P-spline model (TPPS) to account for the local spatial variation exhibited by a response variable measured on a potentially irregular grid of rows and columns of the units. The data may be arranged in sections, for each of which there is a grid and for which the model is to be fitted separately. Also, the rows and columns of a grid are not necessarily one observational unit wide. For TPPS models for which the order of differencing the penalty matrix is two, the an optimal rotation of the null-space eigenvectors of the penalty matrix can be investigated.

No hypothesis testing or comparison of information criteria is made. To use information criteria to decide whether to change the model use chooseSpatialModelOnIC.asrtests.

The model fit supplied in the asrtests.obj should not include terms that will be included in the local spatial model. All spatial model terms are fitted as fixed or random. Consequently, the residual model does not have to be iid.

One or more rows is added for each section to the test.summary data.frame. Convergence and the occurrence of fixed correlations in fitting the model is checked and a note included in the action if there was not. All components of the asrtests.object are updated for the new model.

Usage

## S3 method for class 'asrtests'
addSpatialModel(asrtests.obj, spatial.model = "TPPS", 
                sections = NULL, 
                row.covar = "cRow", col.covar = "cCol", 
                row.factor = "Row", col.factor = "Col", 
                corr.funcs = c("ar1", "ar1"), corr.orders = c(0, 0), 
                row.corrFitfirst = TRUE, 
                allow.corrsJointFit = TRUE, nugget.variance = TRUE, 
                dropFixed = NULL, dropRandom = NULL, 
                nsegs = NULL, nestorder = c(1,1), 
                degree = c(3,3), difforder = c(2,2), 
                usRandLinCoeffs = TRUE, 
                rotateX = FALSE, ngridangles = NULL, 
                which.rotacriterion = "AIC", nrotacores = 1, 
                asreml.option = "grp", tpps4mbf.obj = NULL, 
                allow.unconverged = TRUE, allow.fixedcorrelation = TRUE,
                checkboundaryonly = FALSE, update = TRUE, trace = FALSE, 
                maxit = 30, IClikelihood = "full", which.IC = "AIC", ...)

Arguments

Details

The model to which the spatial models is to be added is supplied in the asrtests.obj. It should not include terms that will be included in the local spatial model. All spatial model terms are fitted as fixed or random. Consequently, the residual model does not have to be iid. The improvement in the fit resulting from the addition of a spatial model to the supplied model is evaluated. Note that the data must be in the order that corresponds to the residual argument with a variable to the right of another variable changes levels in the data frame faster than those of the other variable e.g. Row:Column implies that all levels for Column in consecutive rows of the data.frame with a single Row level.

For the corr spatial model, the default model is an autocorrelation model of order one (ar1) for each dimension. However, any of the single dimension correlation/variance models from asreml can be specified for each dimension, as can no correlation model for a dimension; the models for the two dimensions can differ. Using a forward selection procedure, a series of models are tried, without removing boundary or singular terms, beginning with the addition of row correlation and followed by the addition of column correlation or, if the row.corrFitfirst is set to FALSE, the reverse order. If the fitting of the first-fitted correlation did not result in a model change because the fitting did not converge or correlations were fixed, but the fit of the second correlation was successful, then adding the first correlation will be retried. If one of the metric correlation functions is specified (e.g. exp), then the row.covar or col.covar will be used in the spatial model. However, because the correlations are fitted separately for the two dimensions, the row.factor and col.factor are needed for all models and is used for a dimension that does not involve a correlation/variance function for the fit being performed. Also, the correlation models are fitted as random terms and so the correlation model will include a variance parameter for the grid even when ar1 is used to specify the correlation model, i.e. the model fitted is a variance model and there is no difference between ar1 and ar1v in fitting the model. The variance parameter for this term represents the spatial variance and the fit necessarily includes a nugget term, this being the residual variance. If any correlation is retained in the model, for a section if sections is not NULL, then the need for a nugget term is assessed by fixing the corresponding residual variance to one, unless there are multiple residual variances and these are not related to the sections. Once the fitting of the correlation model has been completed, the rmboundary function will be executed with the checkboundaryonly value supplied in the addSpatialModel.asrtests call. Finally, checking for bound and singular random terms associated with the correlation model and residual terms will be carried out when there are correlation terms in the model and checkboundaryonly has been set to FALSE; as many as possible will be removed from the fitted model, in some cases by fixing variance terms to one.

The tensor-product natural-cubic-smoothing-spline (TPNCSS) spatial model is as described by Verbyla et al. (2018), the tensor-product penalized-cubic-spline (TPPSC2) model with second-order differencing of the penalty is similar to that described by Rodriguez-Alvarez et al. (2018), and the tensor-product, first-difference-penalty, linear spline (TPPSL1) model is amongst those described by Piepho, Boer and Williams (2022). The fixed terms for the spline models are row.covar + col.covar + row.covar:col.covar and the random terms are spl(row.covar) + spl(col.covar) + dev(row.covar) + dev(col.covar) + spl(row.covar):col.covar + row.covar:spl(col.covar) + spl(row.covar):spl(col.covar), except that spl(row.covar) + spl(col.covar) is replaced with spl(row.covar):int(col.covar) + int(row.covar):spl(col.covar) in the TPPSC2 model, where int(.) indicates an intercept or constant value specific to its argument. For TPPSL1 models, the terms spl(row.covar):col.covar + row.covar:spl(col.covar) are omitted, The supplied model should not include any of these terms. However, any fixed or random main-effect Row or Column term that has been included as an initial model for comparison with a spatial model can be removed prior to fitting the spatial model using dropFixed or dropRandom. For the P-spline models with second-order differencing, the model matrices used to fit the pairs of random terms (i) spl(row.covar):int(col.covar) and spl(row.covar):col.covar and (ii) int(row.covar):spl(col.covar) and row.covar:spl(col.covar) are transformed using the spectral decomposition of their penalty matrices. An unstructured variance model is tried for each of these pairs. For TPPSC2, it is also possible to optimize the rotation of the null-space eigenvectors of the penalty matrix for each of these random-term pairs (for more information see Piepho, Boer and Williams, 2022). The optimization is achieved either using an optimizer or takes the form of a search over a grid of rotation angles for a reduced model; the fit of the full model with rotation using the optimal rotation angles will be returned.

The TPPCS and TPP1LS models are fitted using functions from the R package TPSbits authored by Sue Welham (2022). There are two methods for supplying the spline basis information produced by tpsmmb to asreml. The grp method adds it to the data.frame supplied in the data argument of the asreml call. The mbf method creates smaller data.frames with the spline basis information in the same environment as the internal function that calls the spline-fitting function. If it is desired to use in a later session, an asreml function, or asrtests function that calls asreml, (e.g. predict.asreml, predictPlus.asreml, or changeTerms.asrtests) on an asreml.object created using mbf terms, then the mbf data.frames will need to be recreated using makeTPPSplineMats.data.frame in the new session, supplying, if there has been rotation of the penalty matrix eigenvectors, the theta values that are returned as the attribute theta.opt of the asreml.obj.

All models utlize the function changeTerms.asrtests to fit the spatial model. Arguments from tpsmmb and changeTerms.asrtests can be supplied in calls to addSpatialModel.asrtests and will be passed on to the relevant function through the ellipses argument (...).

The data for experiment can be divided sections and the same spatial model fitted separately to each. The fit over all of the sections is assessed. For more detail see sections above.

Each combination of a row.coords and a col.coords does not have to specify a single observation; for example, to fit a local spatial model to the main units of a split-unit design, each combination would correspond to a main unit and all subunits of the main unit would have the same combination.

Value

An asrtests.object containing the components (i) asreml.obj, possibly with attribute theta.opt, (ii) wald.tab, and (iii) test.summary for the model that includes the spatial model, unless the spatial model fails to be fitted when allow.unconverged and/or allow.fixedcorrelation is set to FALSE. If the asrtests.object is the result of fitting a TPPCS model with an exploration of the rotation of the eigenvectors of the penalty matrix for the linear components, then the asreml.obj will have an attribute theta.opt that contains the optimal rotation angles of the eigenvectors.

Author(s)

Chris Brien

References

Piepho, H.-P., Boer, M. P., & Williams, E. R. (2022). Two-dimensional P-spline smoothing for spatial analysis of plant breeding trials. Biometrical Journal, 64, 835-857.

Rodriguez-Alvarez, M. X., Boer, M. P., van Eeuwijk, F. A., & Eilers, P. H. C. (2018). Correcting for spatial heterogeneity in plant breeding experiments with P-splines. Spatial Statistics, 23, 52-71.

Verbyla, A. P., De Faveri, J., Wilkie, J. D., & Lewis, T. (2018). Tensor Cubic Smoothing Splines in Designed Experiments Requiring Residual Modelling. Journal of Agricultural, Biological and Environmental Statistics, 23(4), 478-508.

Welham, S. J. (2022) TPSbits: Creates Structures to Enable Fitting and Examination of 2D Tensor-Product Splines using ASReml-R. Version 1.0.0 https://mmade.org/tpsbits/

See Also

as.asrtests, makeTPPSplineMats.data.frame, addSpatialModelOnIC.asrtests,
chooseSpatialModelOnIC.asrtests, changeModelOnIC.asrtests, changeTerms.asrtests,
rmboundary.asrtests, testranfix.asrtests, testresidual.asrtests, newfit.asreml,
reparamSigDevn.asrtests, changeTerms.asrtests, infoCriteria.asreml

Examples

## Not run: 

data(Wheat.dat)

#Add row and column covariates
Wheat.dat <- within(Wheat.dat, 
                    {
                      cColumn <- dae::as.numfac(Column)
                      cColumn <- cColumn  - mean(unique(cColumn))
                      cRow <- dae::as.numfac(Row)
                      cRow <- cRow - mean(unique(cRow))
                    })

#Fit initial model
current.asr <- asreml(yield ~ Rep + WithinColPairs + Variety, 
                      random = ~ Row + Column,
                      data=Wheat.dat)

#Create an asrtests object, removing boundary terms
current.asrt <- as.asrtests(current.asr, NULL, NULL, 
                            label = "Random Row and Column effects")
current.asrt <- rmboundary(current.asrt)

#Create an asrtests object with a P-spline spatial variation model 
spatial.asrt <- addSpatialModel(current.asrt, spatial.model = "TPPS", 
                                row.covar = "cRow", col.covar = "cColumn",
                                dropRowterm = "Row", dropColterm = "Column",
                                asreml.option = "grp")
infoCriteria(current.asrt$asreml.obj)

#Create an asrtests object with a P-spline spatial variation model 
#that includes rotation of the eigenvectors of the penalty matrix
spatial.asrt <- addSpatialModel(current.asrt, spatial.model = "TPPS", 
                                row.covar = "cRow", col.covar = "cColumn",
                                dropRowterm = "Row", dropColterm = "Column",
                                rotateX = TRUE, 
                                which.rotacriterion = "dev", 
                                nrotacores = parallel::detectCores(), 
                                asreml.option = "mbf")
infoCriteria(current.asrt$asreml.obj)

## End(Not run)

Description

Adds either a correlation, two-dimensional tensor-product natural cubic smoothing spline (TPNCSS), or a two-dimensional tensor-product penalized P-spline model (TPPS) to account for the local spatial variation exhibited by a response variable measured on a potentially irregular grid of rows and columns of the units. The data may be arranged in sections for each of which there is a grid and for which the model is to be fitted separately. Also, the rows and columns of a grid are not necessarily one observational unit wide. The spatial model is only added if the information criterion of the supplied model is decreased with the addition of the local spatial model. For TPPS models for which the order of differencing the penalty matrix is two, the improvement in the fit from rotating the eigenvectors of the penalty matrix can be investigated; if there is no improvement, the unrotated fit will be returned.

A row is added for each section to the test.summary data.frame of the asrtests.object stating whether or not the new model has been swapped for a model in which the spatial model has been add to the supplied model. Convergence and the occurrence of fixed correlations in fitting the model is checked and a note included in the action if there was not. All components of the asrtests.object are updated to exhibit the differences between the supplied and the new model, if a spatial model is added.

Usage

## S3 method for class 'asrtests'
addSpatialModelOnIC(asrtests.obj, spatial.model = "TPPS", 
                    sections = NULL, 
                    row.covar = "cRow", col.covar = "cCol", 
                    row.factor = "Row", col.factor = "Col", 
                    corr.funcs = c("ar1", "ar1"), corr.orders = c(0, 0), 
                    row.corrFitfirst = TRUE, 
                    allow.corrsJointFit = TRUE, nugget.variance = TRUE, 
                    dropFixed = NULL, dropRandom = NULL, 
                    nsegs = NULL, nestorder = c(1,1), 
                    degree = c(3,3), difforder = c(2,2), 
                    usRandLinCoeffs = TRUE, 
                    rotateX = FALSE, ngridangles = NULL, 
                    which.rotacriterion = "AIC", nrotacores = 1, 
                    asreml.option = "grp", tpps4mbf.obj = NULL, 
                    allow.unconverged = TRUE, allow.fixedcorrelation = TRUE,
                    checkboundaryonly = FALSE, update = TRUE, trace = FALSE, 
                    maxit = 30, IClikelihood = "full", which.IC = "AIC", ...)

Arguments

Details

A fitted spatial model is only returned if it improves the fit over and above that of achieved with the model fit supplied in the asrtests.obj. To fit the spatial model without any hypotheses testing or comparison of information criteria use addSpatialModel.asrtests. The model fit supplied in the asrtests.obj should not include terms that will be included in the local spatial model. All spatial model terms are fitted as fixed or random. Consequently, the residual model does not have to be iid. Note that the data must be in the order that corresponds to the residual argument with a variable to the right of another variable changes levels in the data frame faster than those of the other variable e.g. Row:Column implies that all levels for Column in consecutive rows of the data.frame with a single Row level.

For the corr spatial model, the default model is an autocorrelation model of order one (ar1) for each dimension. However, any of the single dimension correlation/variance models from asreml can be specified for each dimension, as can no correlation model for a dimension; the models for the two dimensions can differ. Using a forward selection procedure, a series of models are tried, without removing boundary or singular terms, beginning with the addition of row correlation and followed by the addition of column correlation or, if the row.corrFitfirst is set to FALSE, the reverse order. If the fitting of the first-fitted correlation did not result in a model change because the fitting did not converge or correlations were fixed, but the fit of the second correlation was successful, then adding the first correlation will be retried. If one of the metric correlation functions is specified (e.g. exp), then the row.covar or col.covar will be used in the spatial model. However, because the correlations are fitted separately for the two dimensions, the row.factor and col.factor are needed for all models and is used for a dimension that does not involve a correlation/variance function for the fit being performed. Also, the correlation models are fitted as random terms and so the correlation model will include a variance parameter for the grid even when ar1 is used to specify the correlation model, i.e. the model fitted is a variance model and there is no difference between ar1 and ar1v in fitting the model. The variance parameter for this term represents the spatial variance and the fit necessarily includes a nugget term, this being the residual variance. If any correlation is retained in the model, for a section if sections is not NULL, then the need for a nuggest term is assessed by fixing the corresponding residual variance to one, unless there are multiple residual variances and these are not related to the sections. Once the fitting of the correlation model has been completed, the rmboundary function will be executed with the checkboundaryonly value supplied in the addSpatialModelOnIC.asrtests call. Finally, checking for bound and singular random terms associated with the correlation model and residual terms will be carried out when there are correlation terms in the model and checkboundaryonly has been set to FALSE; as many as possible will be removed from the fitted model, in some cases by fixing variance terms to one.

The tensor-product natural-cubic-smoothing-spline (TPNCSS) spatial model is as described by Verbyla et al. (2018), the tensor-product penalized-cubic-spline (TPPSC2) model with second-order differencing of the penalty is similar to that described by Rodriguez-Alvarez et al. (2018), and the tensor-product, first-difference-penalty, linear spline (TPPSL1) model is amongst those described by Piepho, Boer and Williams (2022). The fixed terms for the spline models are row.covar + col.covar + row.covar:col.covar and the random terms are spl(row.covar) + spl(col.covar) + dev(row.covar) + dev(col.covar) + spl(row.covar):col.covar + row.covar:spl(col.covar) + spl(row.covar):spl(col.covar), except that spl(row.covar) + spl(col.covar) is replaced with spl(row.covar):int(col.covar) + int(row.covar):spl(col.covar) in the TPPSC2 model, where int(.) indicates an intercept or constant value specific to its argument. For TPPSL1 models, the terms spl(row.covar):col.covar + row.covar:spl(col.covar) are omitted, The supplied model should not include any of these terms. However, any fixed or random main-effect Row or Column term that has been included as an initial model for comparison with a spatial model can be removed prior to fitting the spatial model using dropFixed or dropRandom. For the P-spline models with second-order differencing, the model matrices used to fit the pairs of random terms (i) spl(row.covar):int(col.covar) and spl(row.covar):col.covar and (ii) int(row.covar):spl(col.covar) and row.covar:spl(col.covar) are transformed using the spectral decomposition of their penalty matrices. An unstructured variance model is tried for each of these pairs and retained if it improves the fit. For TPPSC2, it is also possible to optimize the rotation of the null-space eigenvectors of the penalty matrix for each of these random-term pairs (for more information see Piepho, Boer and Williams, 2022). The optimization is achieved either using an optimizer or takes the form of a search over a grid of rotation angles for a reduced model; the fit of the full model with rotation using the optimal rotation angles will only be returned if it improves on the fit of the full, unrotated model.

The TPPCS and TPP1LS models are fitted using functions from the R package TPSbits authored by Sue Welham (2022). There are two methods for supplying the spline basis information produced by tpsmmb to asreml. The grp method adds it to the data.frame supplied in the data argument of the asreml call. The mbf method creates smaller data.frames with the spline basis information in the same environment as the internal function that calls the spline-fitting function. If it is desired to use in a later session, an asreml function, or asrtests function that calls asreml, (e.g. predict.asreml, predictPlus.asreml, or changeTerms.asrtests) on an asreml.object created using mbf terms, then the mbf data.frames will need to be recreated using makeTPPSplineMats.data.frame in the new session, supplying, if there has been rotation of the penalty matrix eigenvectors, the theta values that are returned as the attribute theta.opt of the asreml.obj.

All models utlize the function changeModelOnIC.asrtests to assess the model fit, the information criteria used in assessing the fit being calculated using infoCriteria. Any bound terms are removed from the model. Arguments from tpsmmb and changeModelOnIC.asrtests can be supplied in calls to addSpatialModelOnIC.asrtests and will be passed on to the relevant function through the ellipses argument (...).

The data for experiment can be divided sections and the same spatial model fitted separately to each. The fit over all of the sections is assessed. For more detail see sections above.

Each combination of a row.coords and a col.coords does not have to specify a single observation; for example, to fit a local spatial model to the main units of a split-unit design, each combination would correspond to a main unit and all subunits of the main unit would have the same combination.

Value

An asrtests.object containing the components (i) asreml.obj, possibly with attribute theta.opt, (ii) wald.tab, and (iii) test.summary for the model whose fit has the smallest information criterion between the supplied and spatial model. The values of the degrees of freedom and the information criteria in the test.summary are differences between those of the changed model and those of the model supplied to addSpatialModelOnIC. If the asrtests.object is the result of fitting a TPPCS model with an exploration of the rotation of the eigenvectors of the penalty matrix for the linear components, then the asreml.obj will have an attribute theta.opt that contains the optimal rotation angles of the eigenvectors.

Author(s)

Chris Brien

References

Piepho, H.-P., Boer, M. P., & Williams, E. R. (2022). Two-dimensional P-spline smoothing for spatial analysis of plant breeding trials. Biometrical Journal, 64, 835-857.

Rodriguez-Alvarez, M. X., Boer, M. P., van Eeuwijk, F. A., & Eilers, P. H. C. (2018). Correcting for spatial heterogeneity in plant breeding experiments with P-splines. Spatial Statistics, 23, 52-71.

Verbyla, A. P., De Faveri, J., Wilkie, J. D., & Lewis, T. (2018). Tensor Cubic Smoothing Splines in Designed Experiments Requiring Residual Modelling. Journal of Agricultural, Biological and Environmental Statistics, 23(4), 478-508.

Welham, S. J. (2022) TPSbits: Creates Structures to Enable Fitting and Examination of 2D Tensor-Product Splines using ASReml-R. Version 1.0.0 https://mmade.org/tpsbits/

See Also

as.asrtests, makeTPPSplineMats.data.frame, addSpatialModel.asrtests,
chooseSpatialModelOnIC.asrtests, changeModelOnIC.asrtests, changeTerms.asrtests,
rmboundary.asrtests, testranfix.asrtests, testresidual.asrtests, newfit.asreml,
reparamSigDevn.asrtests, changeTerms.asrtests, infoCriteria.asreml

Examples

## Not run: 

data(Wheat.dat)

#Add row and column covariates
Wheat.dat <- within(Wheat.dat, 
                    {
                      cColumn <- dae::as.numfac(Column)
                      cColumn <- cColumn  - mean(unique(cColumn))
                      cRow <- dae::as.numfac(Row)
                      cRow <- cRow - mean(unique(cRow))
                    })

#Fit initial model
current.asr <- asreml(yield ~ Rep + WithinColPairs + Variety, 
                      random = ~ Row + Column,
                      data=Wheat.dat)

#Create an asrtests object, removing boundary terms
current.asrt <- as.asrtests(current.asr, NULL, NULL, 
                            label = "Random Row and Column effects")
current.asrt <- rmboundary(current.asrt)

current.asrt <- addSpatialModelOnIC(current.asrt, spatial.model = "TPPS", 
                                    row.covar = "cRow", col.covar = "cColumn",
                                    dropRowterm = "Row", dropColterm = "Column",
                                    asreml.option = "grp")
infoCriteria(current.asrt$asreml.obj)

## End(Not run)

Description

A row that summarizes the result of a proposed change to a model is added to a test.summary data.frame. Only the values of those arguments for which there are columns in test.summary will be included in the row.

Usage

addto.test.summary(test.summary, terms, DF = 1, denDF = NA, 
                          p = NA, AIC = NA, BIC = NA,
                          action = "Boundary")

Arguments

Value

A data.frame.

Author(s)

Chris Brien

See Also

asremlPlus-package, asrtests.object, print.test.summary

Examples

## Not run: 
data(Wheat.dat)

## Fit an autocorrelation model
ar1.asr <- asreml(yield ~ Rep + WithinColPairs + Variety, 
                  random = ~ Row + Column + units,
                  residual = ~ ar1(Row):ar1(Column), 
                  data=Wheat.dat)
ar1.asrt <- as.asrtests(ar1.asr, NULL, NULL, 
                        label = "Autocorrelation model")
ar1.asrt <- rmboundary.asrtests(ar1.asrt)

## Fit a tensor spline
Wheat.dat <- within(Wheat.dat, 
                    {
                      cRow <- dae::as.numfac(Row)
                      cRow <- cRow - mean(unique(cRow))
                      cColumn <- dae::as.numfac(Column)
                      cColumn <- cColumn - mean(unique(cColumn))
                    })
ts.asr <- asreml(yield ~ Rep + cRow + cColumn + WithinColPairs + 
                         Variety, 
                  random = ~ spl(cRow) + spl(cColumn) + 
                             dev(cRow) + dev(cColumn) + 
                             spl(cRow):cColumn + cRow:spl(cColumn) + 
                             spl(cRow):spl(cColumn),
                  residual = ~ Row:Column, 
                  data=Wheat.dat)
ts.asrt <- as.asrtests(ts.asr, NULL, NULL, 
                       label = "Tensor spline model")
ts.asrt <- rmboundary.asrtests(ts.asrt)
ar1.ic <- infoCriteria(ar1.asrt$asreml.obj)
ts.ic <- infoCriteria(ts.asrt$asreml.obj)
if (ar1.ic$AIC < ts.ic$AIC)
{
  ic.diff <- ar1.ic - ts.ic
  new.asrt <- ar1.asrt 
  new.asrt$test.summary <- addto.test.summary(ar1.asrt$test.summary, 
                                              terms = "Compare ar1 to ts", 
                                              DF = ic.diff$varDF, 
                                              AIC = ic.diff$AIC, BIC = ic.diff$BIC, 
                                              action = "Chose ar1")
} else
{
  ic.diff <- ts.ic - ar1.ic
  new.asrt <- ts.asrt
  new.asrt$test.summary <- addto.test.summary(ts.asrt$test.summary, 
                                              terms = "Compare ar1 to ts", 
                                              DF = ic.diff$varDF, 
                                              AIC = ic.diff$AIC, BIC = ic.diff$BIC, 
                                              action = "Chose ts")
}

## End(Not run)

Description

Uses supplied predictions and standard errors of pairwise differences, or the variance matrix of predictions to form, in an alldiffs.object, for those components not already present, (i) a table of all pairwise differences of the predictions, (ii) the p-value of each pairwise difference, and (iii) the minimum, mean, maximum and accuracy of LSD values. Predictions that are aliased (or inestimable) are removed from the predictions component of the alldiffs.object and standard errors of differences involving them are removed from the sed component.

If necessary, the order of the columns of the variables in the predictions component are changed to be the initial columns of the predictions.frame and to match their order in the classify. Also, the rows of predictions component are ordered so that they are in standard order for the variables in the classify. That is, the values of the last variable change with every row, those of the second-last variable only change after all the values of the last variable have been traversed; in general, the values of a variable are the same for all the combinations of the values to the variables to its right in the classify. The sortFactor or sortOrder arguments can be used to order of the values for the classify variables, which is achieved using sort.alldiffs.

Each p-value is computed as the probability of a t-statistic as large as or larger than the absolute value of the observed difference divided by its standard error. The p-values are stored in the p.differences component. The degrees of freedom of the t-distribution is the degrees of freedom stored in the tdf attribute of the alldiffs.object. This t-distribution is also used in calculating the LSD statistics stored in the LSD component of the alldiffs.object.

Usage

## S3 method for class 'data.frame'
allDifferences(predictions, classify, vcov = NULL, 
               differences = NULL, p.differences = NULL, sed = NULL, 
               LSD = NULL, LSDtype = "overall", LSDsupplied = NULL, 
               LSDby = NULL, LSDstatistic = "mean", 
               LSDaccuracy = "maxAbsDeviation", 
               retain.zeroLSDs = FALSE, 
               zero.tolerance = .Machine$double.eps ^ 0.5, 
               backtransforms = NULL, 
               response = NULL, response.title = NULL, 
               term = NULL, tdf = NULL,  
               x.num = NULL, x.fac = NULL,  
               level.length = NA, 
               pairwise = TRUE, alpha = 0.05,
               transform.power = 1, offset = 0, scale = 1, 
               transform.function = "identity", 
               inestimable.rm = TRUE,
               sortFactor = NULL, sortParallelToCombo = NULL, 
               sortNestingFactor = NULL, sortOrder = NULL, 
               decreasing = FALSE, ...)

Arguments

Value

An alldiffs.object with components predictions, vcov, differences, p.differences sed, and LSD.

The name of the response, the response.title, the term, the classify, tdf, alpha, sortFactor and the sortOrder will be set as attributes to the object. Note that the classify in an alldiffs.object is based on the variables indexing the predictions, which may differ from the classify used to obtain the original predictions (for example, when the alldiffs.objects stores a linear transformation of predictions.

Also, see predictPlus.asreml for more information.

Author(s)

Chris Brien

See Also

asremlPlus-package, as.alldiffs, as.predictions.frame, sort.alldiffs, subset.alldiffs,
print.alldiffs, renewClassify.alldiffs, redoErrorIntervals.alldiffs,
recalcLSD.alldiffs, pickLSDstatistics.alldiffs, plotPredictions.data.frame,
predictPlus.asreml, predictPresent.asreml

Examples

data(Oats.dat)
  
  ## Use asreml to get predictions and associated statistics

  ## Not run: 
  m1.asr <- asreml(Yield ~ Nitrogen*Variety, 
                   random=~Blocks/Wplots,
                   data=Oats.dat)
  current.asrt <- as.asrtests(m1.asr)
  Var.pred <- asreml::predict.asreml(m1.asr, classify="Nitrogen:Variety", 
                                      sed=TRUE)
  if (getASRemlVersionLoaded(nchar = 1) == "3")
    Var.pred <- Var.pred$predictions
  Var.preds <- Var.pred$pvals
  Var.sed <- Var.pred$sed
  Var.vcov <- NULL
  wald.tab <-  current.asrt$wald.tab
  den.df <- wald.tab[match("Variety", rownames(wald.tab)), "denDF"]
  
## End(Not run)

  ## Use lmerTest and emmmeans to get predictions and associated statistics
  if (requireNamespace("lmerTest", quietly = TRUE) & 
      requireNamespace("emmeans", quietly = TRUE))
  {
    m1.lmer <- lmerTest::lmer(Yield ~ Nitrogen*Variety + (1|Blocks/Wplots),
                              data=Oats.dat)
    Var.emm <- emmeans::emmeans(m1.lmer, specs = ~ Nitrogen:Variety)
    Var.preds <- summary(Var.emm)
    den.df <- min(Var.preds$df)
    ## Modify Var.preds to be compatible with a predictions.frame
    Var.preds <- as.predictions.frame(Var.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
    Var.vcov <- vcov(Var.emm)
    Var.sed <- NULL
  }

  ## Use the predictions obtained with either asreml or lmerTest
  if (exists("Var.preds"))
  {
    ## Order the Varieties in decreasing order for the predictions values in the 
    ## first N level 
    Var.diffs <- allDifferences(predictions = Var.preds, 
                                classify = "Nitrogen:Variety", 
                                sed = Var.sed, vcov = Var.vcov, tdf = den.df,
                                sortFactor = "Variety", decreasing = TRUE)
    print.alldiffs(Var.diffs, which="differences")
  
    ## Change the order of the factors in the alldiffs object and reorder components
    Var.reord.diffs <- allDifferences(predictions = Var.preds,
                                classify = "Variety:Nitrogen", 
                                sed = Var.sed, vcov = Var.vcov, tdf = den.df)
    print.alldiffs(Var.reord.diffs, which="predictions")
  }

Description

An object of S3-class alldiffs that stores the predictions for a model, along with supplied statistics for all pairwise differences. While alldiffs.object can be constructed by defining a list with the appropriate components, it can be formed by passing the components to as.alldiffs, or from a predictions data.frame using allDifferences.data.frame.

as.alldiffs is function that assembles an object of this class from supplied components.

is.alldiffs is the membership function for this class; it tests that an object is of class alldiffs.

validAlldiffs(object) can be used to test the validity of an object with this class.

allDifferences.data.frame is the function that constructs an object of this class by calculating components from statistics supplied via its arguments and then using as.alldiffs to make the object.

Value

A list of class alldiffs containing the following components: predictions, vcov, differences, p.differences, sed, LSD and backtransforms. Except for predictions, the components are optional and can be set to NULL.

An alldiffs.object also has attributes response, response.title, term, classify, tdf, alpha, sortFactor and sortOrder, which may be set to NULL.

The details of the components are as follows:

1. predictions: A predictions.frame, being a data.frame beginning with the variables classifying the predictions, in the same order as in the classify, and also containing columns named predicted.value, standard.error and est.status; each row contains a single predicted value. The number of rows should equal the number of unique combinations of the classify variables and will be in standard order for the classify variables. That is, the values of the last variable change with every row, those of the second-last variable only change after all the values of the last variable have been traversed; in general, the values of a variable are the same for all the combinations of the values to the variables to its right in the classify.

   The data.frame may also include columns for the lower and upper values of error intervals, either standard error, confidence or half-LSD intervals. The names of these columns will consist of three parts separated by full stops: 1) the first part will be lower or upper; 2) the second part will be one of Confidence, StandardError or halfLeastSignificant; 3) the third component will be limits.

   Note that the names standard.error and est.status have been changed to std.error and status in the pvals component produced by asreml-R4; if the new names are in the data.frame supplied to predictions, they will be returned to the previous names.

2. differences: A matrix containing all pairwise differences between the predictions; it should have the same number of rows and columns as there are rows in predictions.

3. p.differences: A matrix containing p-values for all pairwise differences between the predictions; each p-value is computed as the probability of a t-statistic as large as or larger than the observed difference divided by its standard error. The degrees of freedom of the t distribution for computing it are computed as the denominator degrees of freedom of the F value for the fixed term, if available; otherwise, the degrees of freedom stored in the attribute tdf are used; the matrix should be of the same size as that for differences.

4. sed: A matrix containing the standard errors of all pairwise differences between the predictions; they are used in computing the p-values in p.differences.

5. vcov: A matrix containing the variance matrix of the predictions; it is used in computing the variance of linear transformations of the predictions.

6. LSD: An LSD.frame containing (i) c, the number of pairwise predictions comparisons for each LSD value and the mean, minimum, maximum and assigned LSD, (ii) the column accuracyLSD that gives a measure of the accuracy of the assigned LSD. given the variation in LSD values, and (iii) the columns false.pos and false.neg that contain the number of false positives and negatives if the assignedLSD value(s) is(are) used to determine the significance of the pairwise predictions differences. The LSD values in the assignedLSD column is used to determine the significance of pairwise differences that involve predictions for the combination of levels given by a row name. The value in the assignedLSD column is specified using the LSDstatistic argument.

7. backtransforms: When the response values have been transformed for analysis, a data.frame containing the backtransformed values of the predicted values is added to the alldiffs.object. This data.frame is consistent with the predictions component, except that the column named predicted.value is replaced by one called backtransformed.predictions. Any error.interval values will also be the backtransformed values. Each row contains a single predicted value.

The details of the attributes of an alldiffs.object are:

1. response: A character specifying the response variable for the predictions.

2. response.title: A character specifying the title for the response variable for the predictions.

3. term: A character giving the variables that define the term that was fitted using asreml and that corresponds to classify. It is often the same as classify.

4. classify: A character giving the variables that define the margins of the multiway table used in the prediction. Multiway tables are specified by forming an interaction type term from the classifying variables, that is, separating the variable names with the : operator.

5. tdf: An integer specifying the degrees of freedom of the standard error. It is used as the degrees of freedom for the t-distribution on which p-values and confidence intervals are based.

6. alpha: An integer specifying the significance level. It is used as the significance level calculating LSDs.

7. LSDtype: If the LSD component is not NULL then LSDtype is added as an attribute. A character nominating the type of grouping of seds to be used in combining LSDs.

8. LSDby: If the LSD component is not NULL then LSDby is added as an attribute. A character vector containing the names of the factors and numerics within whose combinations the LSDs are to be summarized.

9. LSDstatistic: If the LSD component is not NULL then LSDstatistic is added as an attribute. A character nominating what statistic to use in summarizing a set of LSDs.

10. LSDaccuracy: If the LSD component is not NULL then LSDaccuracy is added as an attribute. A character nominating the method of calculating a measure of the accuracy of the LSDs stored in the assignedLSD column of the LSD.frame.

11. sortFactor: factor that indexes the set of predicted values that determined the sorting of the components.

12. sortOrder: A character vector that is the same length as the number of levels for sortFactor in the predictions component of the alldiffs.object. It specifies the order of the levels in the reordered components of the alldiffs.object.

The following creates a sortOrder vector levs for factor f based on the values in x:
levs <- levels(f)[order(x)].

See predictPlus.asreml for more information.

Author(s)

Chris Brien

See Also

is.alldiffs, as.alldiffs, validAlldiffs, allDifferences.data.frame

Examples

data(Oats.dat)

  ## Use asreml to get predictions and associated statistics

  ## Not run: 
  m1.asr <- asreml(Yield ~ Nitrogen*Variety, 
                   random=~Blocks/Wplots,
                   data=Oats.dat)
  current.asrt <- as.asrtests(m1.asr)
  Var.pred <- asreml::predict.asreml(m1.asr, classify="Nitrogen:Variety", 
                                      sed=TRUE)
  if (getASRemlVersionLoaded(nchar = 1) == "3")
    Var.pred <- Var.pred$predictions
  Var.preds <- Var.pred$pvals
  Var.sed <- Var.pred$sed
  Var.vcov <- NULL
  
## End(Not run)
  
  ## Use lmerTest and emmmeans to get predictions and associated statistics
  if (requireNamespace("lmerTest", quietly = TRUE) & 
      requireNamespace("emmeans", quietly = TRUE))
  {
    m1.lmer <- lmerTest::lmer(Yield ~ Nitrogen*Variety + (1|Blocks/Wplots),
                              data=Oats.dat)
    Var.emm <- emmeans::emmeans(m1.lmer, specs = ~ Nitrogen:Variety)
    Var.preds <- summary(Var.emm)
    den.df <- min(Var.preds$df)
    ## Modify Var.preds to be compatible with a predictions.frame
    Var.preds <- as.predictions.frame(Var.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
    Var.vcov <- vcov(Var.emm)
    Var.sed <- NULL
  }

  ## Use the predictions obtained with either asreml or lmerTest
  if (exists("Var.preds"))
  {
    ## Form an all.diffs object
     Var.diffs <- as.alldiffs(predictions = Var.preds, classify = "Nitrogen:Variety", 
                              sed = Var.sed, vcov = Var.vcov, tdf = den.df)

    ## Check the class and validity of the alldiffs object
    is.alldiffs(Var.diffs)
    validAlldiffs(Var.diffs)
  }

Description

Applies the angular transformation to numeric values. It is given by $\sin^{-1}(\sqrt{proportions})$

Usage

angular(proportions, n)

Arguments

Value

A numeric.

Author(s)

Chris Brien

See Also

angular.mod, powerTransform.

Examples

n <-25
y <- rbinom(10, n, 0.5)
y <- c(y,0,n)
p <- y/n
p.ang <- angular(p, n)

Description

Applies the angular transformation to a vector of counts. A modified transformation is used that is appropriate when N < 50 and the proportion is not between 0.3 and 0.7. The transformation is given by $\sin^{-1}{\frac{count + 0.375}{n + 0.75}}$.

Usage

angular.mod(count, n)

Arguments

Value

A numeric vector.

Author(s)

Chris Brien

See Also

angular, powerTransform.

Examples

n <-25
y <- rbinom(10, n, 0.5)
y <- c(y,0,n)
p.ang.mod <- angular.mod(y, n)

Forms an alldiffs.object from the supplied predictions, along with those statistics, associated with the predictions and their pairwise differences, that have been supplied.

Description

Creates an alldiffs.object that consists of a list containing the following components: predictions, vcov, differences, p.differences, sed, LSD and backtransforms. Predictions must be supplied to the function while the others will be set only if they are supplied; those not supplied are set to NULL. It also has attributes response, response.title, term, classify, tdf, tdf, alpha, sortFactor and sortOrder. which will be set to the values supplied or NULL if none are supplied.

Usage

as.alldiffs(predictions, vcov = NULL, differences = NULL, 
            p.differences = NULL, sed = NULL, LSD = NULL, 
            backtransforms = NULL, 
            response = NULL, response.title = NULL, 
            term = NULL, classify = NULL, 
            tdf = NULL, alpha = 0.05, 
            sortFactor = NULL, sortOrder = NULL)

Arguments

Value

An S3-class alldiffs.object. Also, see predictPlus.asreml for more information.

Author(s)

Chris Brien

See Also

asremlPlus-package, alldiffs.object, is.alldiffs, as.alldiffs, print.alldiffs,
sort.alldiffs, subset.alldiffs, allDifferences.data.frame,
renewClassify.alldiffs, redoErrorIntervals.alldiffs, recalcLSD.alldiffs,
predictPlus.asreml, plotPredictions.data.frame, predictPresent.asreml

Examples

data(Oats.dat)
  
  ## Use asreml to get predictions and associated statistics

  ## Not run: 
  m1.asr <- asreml(Yield ~ Nitrogen*Variety, 
                   random=~Blocks/Wplots,
                   data=Oats.dat)
  current.asrt <- as.asrtests(m1.asr)
  Var.pred <- asreml::predict.asreml(m1.asr, classify="Nitrogen:Variety", 
                                      sed=TRUE)
  if (getASRemlVersionLoaded(nchar = 1) == "3")
    Var.pred <- Var.pred$predictions
  Var.preds <- Var.pred$pvals
  Var.sed <- Var.pred$sed
  Var.vcov <- NULL
  
## End(Not run)
  
  ## Use lmerTest and emmmeans to get predictions and associated statistics
  if (requireNamespace("lmerTest", quietly = TRUE) & 
      requireNamespace("emmeans", quietly = TRUE))
  {
    m1.lmer <- lmerTest::lmer(Yield ~ Nitrogen*Variety + (1|Blocks/Wplots),
                              data=Oats.dat)
    Var.emm <- emmeans::emmeans(m1.lmer, specs = ~ Nitrogen:Variety)
    Var.preds <- summary(Var.emm)
    den.df <- min(Var.preds$df)
    ## Modify Var.preds to be compatible with a predictions.frame
    Var.preds <- as.predictions.frame(Var.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
    Var.vcov <- vcov(Var.emm)
    Var.sed <- NULL
  }

  ## Use the predictions obtained with either asreml or lmerTest
  if (exists("Var.preds"))
  {
    ## Form an all.diffs object
     Var.diffs <- as.alldiffs(predictions = Var.preds, classify = "Nitrogen:Variety", 
                              sed = Var.sed, vcov = Var.vcov, tdf = den.df)

    ## Check the class and validity of the alldiffs object
    is.alldiffs(Var.diffs)
    validAlldiffs(Var.diffs)
  }

Description

An asrtests.object that is a list consisting of the components asreml.obj, wald.tab and test.summary.

A call to as.asrtests with test.summary = NULL re-initializes the test.summary
data.frame.

If there is no wald.tab, wald.asreml is called. In all cases, recalcWaldTab is called and any changes made as specified by the recalcWaldTab arguments supplied via ....

The label argument can be used to include an entry in test.summary for the starting model. If a label is included, (i) the information criteria calculated using the asreml.obj will be added to the test.summary, if IClikelihood is not set to none and (ii) the number of variance parameters is included in the denDF column, if IClikelihood is set to none.

Usage

as.asrtests(asreml.obj, wald.tab = NULL, test.summary = NULL, 
            denDF = "numeric", label = NULL, 
            IClikelihood = "none", bound.exclusions = c("F","B","S","C"), ...)

Arguments

Value

An object of S3-class asrtests.

Author(s)

Chris Brien

References

Kenward, M. G., & Roger, J. H. (1997). Small sample inference for fixed effects from restricted maximum likelihood. Biometrics, 53, 983-997.

See Also

asremlPlus-package, is.alldiffs, as.alldiffs, recalcWaldTab,
testranfix.asrtests, chooseModel.asrtests, rmboundary.asrtests,
reparamSigDevn.asrtests

Examples

## Not run: 
data(Wheat.dat)

# Fit initial model
current.asr <- asreml(yield ~ Rep + WithinColPairs + Variety, 
                      random = ~ Row + Column + units,
                      residual = ~ ar1(Row):ar1(Column), 
                      data=Wheat.dat)

# Load current fit into an asrtests object
current.asrt <- as.asrtests(current.asr, NULL, NULL)

# Check for and remove any boundary terms
current.asrt <- rmboundary(current.asrt)

## End(Not run)

Forms a predictions.frame from a data.frame, ensuring that the correct columns are present.

Description

Creates a predictions.frame from a data.frame by adding the class predictions.frame to it, and renaming the columns containing the predictions, se, est.status and error.intervals.

Usage

as.predictions.frame(data, classify = NULL, 
                     predictions = NULL, se = NULL, est.status = NULL, 
                     interval.type = NULL, interval.names = NULL)

Arguments

Value

An S3-class predictions.frame.

Author(s)

Chris Brien

See Also

asremlPlus-package, predictions.frame, is.predictions.frame, predictions.frame,
validPredictionsFrame

Examples

data(Oats.dat)
  
  ## Use asreml to get predictions and associated statistics

  ## Not run: 
  m1.asr <- asreml(Yield ~ Nitrogen*Variety, 
                   random=~Blocks/Wplots,
                   data=Oats.dat)
  current.asrt <- as.asrtests(m1.asr)
  Var.pred <- asreml::predict.asreml(m1.asr, classify="Nitrogen:Variety", 
                                      sed=TRUE)
  if (getASRemlVersionLoaded(nchar = 1) == "3")
    Var.pred <- Var.pred$predictions
 #Form predictions.frame changing asreml-R4 names to the standard names, if these are present
 Var.preds <- as.predictions.frame(Var.pred$pvals, se = "std.error", 
                                    est.status = "status")
  
## End(Not run)
  
  ## Use lmerTest and emmmeans to get predictions and associated statistics
  if (requireNamespace("lmerTest", quietly = TRUE) & 
      requireNamespace("emmeans", quietly = TRUE))
  {
    m1.lmer <- lmerTest::lmer(Yield ~ Nitrogen*Variety + (1|Blocks/Wplots),
                              data=Oats.dat)
    Var.emm <- emmeans::emmeans(m1.lmer, specs = ~ Nitrogen:Variety)
    Var.preds <- summary(Var.emm)
    Var.preds <- as.predictions.frame(Var.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
  }

  ## Check the class and validity of the alldiffs object
  if (exists("Var.preds"))
  {
    is.predictions.frame(Var.preds)
    validPredictionsFrame(Var.preds)
  }

Description

These functions have been renamed and deprecated in asremlPlus:

1. addrm.terms.asreml and addrm.terms.asrtests -> changeTerms.asrtests,

2. alldiffs -> as.alldiffs,

3. asrtests-> as.asrtests,

4. choose.model.asreml and choose.model.asrtests -> chooseModel.asrtests,

5. facRecode and facRecode.alldiffs -> facRecast.alldiffs,

6. info.crit and info.crit.asreml -> infoCriteria.asreml,

7. newrcov.asrtests -> changeTerms.asrtests,

8. plotvariofaces.asreml -> plotVariofaces.data.frame,

9. power.transform -> powerTransform,

10. predictiondiffs.asreml -> allDifferences.data.frame,

11. predictionplot.asreml -> plotPredictions.data.frame,

12. predictparallel.asreml -> predictPlus.asreml,

13. pred.present.asreml -> predictPresent.asreml,

14. recalc.wald.tab.asreml and recalc.wald.tab.asrtests -> recalcWaldTab.asrtests,

15. reorderClassify and reorderClassify.alldiffs -> renewClassify.alldiffs,

16. reml.lrt and reml.lrt.asreml -> REMLRT.asreml,

17. rmboundary.asreml -> rmboundary.asrtests,

18. setvarianceterms.asreml -> setvarianceterms.call,

19. sig.devn.reparam.asreml and sig.devn.reparam.asrtests -> reparamSigDevn.asrtests,

20. testranfix.asreml -> testranfix.asrtests,

21. testrcov.asreml and testrcov.asrtests -> testresidual.asrtests,

22. testswapran.asreml -> testswapran.asrtests

Usage

addrm.terms.asreml(...)
addrm.terms.asrtests(...)
alldiffs(...)
asrtests(...)
choose.model.asreml(...)
choose.model.asrtests(...)
facRecode(...)
facRecode.alldiffs(...)
info.crit(...)
info.crit.asreml(...)
newrcov.asrtests(...)
plotvariofaces.asreml(...)
power.transform(...)
predictiondiffs.asreml(...)
predictionplot.asreml(...)
predictparallel.asreml(...)
pred.present.asreml(...)
recalc.wald.tab.asreml(...)
recalc.wald.tab.asrtests(...)
reml.lrt(...)
reml.lrt.asreml(...)
## S3 method for class 'alldiffs'
reorderClassify(...)
## S3 method for class 'asreml'
rmboundary(...)
setvarianceterms.asreml(...)
sig.devn.reparam.asreml(...)
sig.devn.reparam.asrtests(...)
testranfix.asreml(...)
testrcov.asreml(...)
testrcov.asrtests(...)
## S3 method for class 'asreml'
testswapran(...)

Arguments

Author(s)

Chris Brien

Description

The intermittent, randomly-presented, startup tips.

Startup tips

Need help? The manual is a vignette and is in the vignettes subdirectory of the package's install directory.

Find out what has changed in asremlPlus: enter news(package = 'asremlPlus').

Need help getting started? Enter vignette(package = 'asremlPlus').

To avoid start-up message that ASReml-R is needed, load asreml before asremlPlus.

The methods for alldiffs and data.frame do not require asreml

Use suppressPackageStartupMessages() to eliminate all package startup messages.

To see all the intermittent, randomly-presented, startup tips enter ?asremlPlusTips.

To install the latest version: go to http://chris.brien.name/rpackages.

For versions between CRAN releases (and more) go to http://chris.brien.name/rpackages.

Author(s)

Chris Brien

Description

An object of S3-class asrtests that contains information derived from the fits of a mixed model using asreml.

as.asrtests is function that makes an object of this class.

is.list is the membership function for this class; it tests that an object is of class list.

validAsrtests can be used to test the validity of an asrtests.object.

Value

A list that contains three components:

1. asreml.obj: an object of class asreml that contains the fit of a model;

2. wald.tab: A data.frame containing a pseudo-anova table for the fixed terms produced by wald.asreml. It has rownames that correspond to the fixed terms that were fitted and four columns. If denominator degrees of freedom were calculated then the columns are DF, denDF, F.inc, Pr; otherwise the columns are Df, Sum of Sq, Wald statistic, and Pr(Chisq).

3. test.summary: A data.frame with columns terms, DF, denDF, p, AIC, BIC and action, each row of which summarizes the results of proposed changes to the fitted model.

   Possible codes for action are: Dropped, Retained, Swapped, Unswapped, Unchanged, Significant, Nonsignificant, Absent, Added, Removed and Boundary. If the either of the models did not converge, unconverged will be added to the code. Unchanged is used when allow.unconverged is FALSE. Note that the logical asreml.obj$converge also reflects whether there is convergence.

   A row is added to the test.summary for each term that is dropped, added or tested or a note that several terms have been added or removed. When values for the AIC and BIC are included in the row, then the DF are the number of fixed parameters in the model and denDF are the numbers of variance parameters. When changeModelOnIC adds a row then the values of the degrees of freedom and information criteria are differences between those for the model that is supplied and the model changed by changeModelOnIC.

Author(s)

Chris Brien

See Also

as.asrtests, as.asrtests, validAsrtests

Description

Extracts the REML log likelihood for two asreml objects and forms the observed REML ratio statistic. It assumes that the second asreml object is the result of fitting a model that is a reduced version of the model for the first object and is considered to the null model. Using the mean and V, nboot bootstrap samples of simulated response values are generated in parallel; that is, ncores cores are used and each is used to generate and analyse a sample. The full and reduced models are fitted to the data and if either analysis fails to converge another sample is generated and analysed using the current core, with a maximum of max.retries attempts to obtain a sample that converges for both analysis. Thus the maximum number of data sets that will be generated is nboot * max.retries. If a bootstrap sample converges for both analyses, the REML ratio test statistic is formed for it. The p-value is then calculated as $(k + 1) / (b + 1)$ where $k$ is the number of simulated ratio test statistics greater than the observed test statistic and $s$ is the number of bootstrap samples that were returned.

The function checks that the models do not differ in either their fixed or sparse models. It also check the difference in the number of variance parameters between the two fits to the models, taking into account the bound.exclusions.

Usage

## S3 method for class 'asreml'
bootREMLRT(h0.asreml.obj, h1.asreml.obj, 
           nboot = 100, max.retries = 5, seed = NULL, 
           means=NULL, V = NULL, extra.matrix = NULL, ignore.terms = NULL, 
           fixed.spline.terms = NULL, 
           bound.exclusions = c("F","B","S","C"), 
           tolerance = 1E-10, update = TRUE, trace = FALSE, 
           ncores = 2, ...)

Arguments

Value

A list with the following components:

1. REMLRT: the observed REML ratio statistic.

2. p: the bootstrap p-value for the observed test statistic.

3. DF: the calculated difference in DF for the variance parameters in the two models.

4. totalunconverged: the total number of unconverged analyses over the simulations.

5. REMLRT.sim: a numeric containing the values of the ratio statistics for the simulated data. It has an attribute called na.action that can be retrieved using attr(REMLRT.sim, which = "na.action"); it contains a list of the simulation numbers that were abandoned because max.retries failed to converge for both models.

6. nunconverged: the number of unconverged analyses for each bootstrap sample, the maximum being max.retries.

Note

A bootstrap sample is generated using a multivariate normal distribution with expected value as specified by means and variance matrix given by V. Each simulated sample is analysed according to the reduced model and, provided this analysis converges, according to the full.model. If one of these analyses fails to converge, it is abandoned and another sample is generated for this simulation. As many as max.retries attempts are made to generate a data set for which both analyses converge. If data set that converges for both analyses is not generated for a simulation, NA is returned for that bootstrap sample. Hence, the maximum number of data sets that will be generated is nboot * max.retries and less than nboot samples will be generated if a data set that converges for both analyses is not obtained within max.retries attempts.

If a bootstrap sample converges for both analyses, the REML ratio test statistic is calculated as $2(log(REML)_F - log(REML)_R)$.

The DF is calculated from the information in full.asreml.obj and reduced.asreml.obj. The degrees of freedom are computed as the difference between the two models in the number of variance parameters whose estimates do not have a code for bound specified in bound.exclusions.

If ASReml-R version 4 is being used then the codes specified in bound.exclusions are not restricted to a subset of the default codes, but a warning is issued if a code other than these is specified. For ASReml-R version 3, only a subset of the default codes are allowed: F (Fixed), B (Boundary), C (Constrained) and S (Singular).

Author(s)

Chris Brien

See Also

REMLRT.asreml, infoCriteria.asreml, newfit.asreml, testranfix.asrtests

Examples

## Not run: 
    bootREMLRT(ICV.max, ICV.red, ncores = parallel::detectCores())

## End(Not run)

Description

Uses information criteria to decide whether to change the fitted model stored in the supplied asrtests.object according to the specified modifications. The function changeTerms is used to change the model. Thus, the model can be modified using a combination of adding and removing sets of terms from one or both of the fixed or random models, replacing the residual model and changing the bounds and/or initial values of some terms. The model will be unchanged if terms specified in dropFixed or dropRandom are not in the fitted model.

A row is added to the test.summary data.frame of the asrtests.object using the supplied label and stating whether or not the new model has been swapped for the supplied model. Convergence in fitting the model is checked and a note included in the action if there was not. All components of the asrtests.object are updated to exhibit the differences between the supplied and new models.

To obtain a list of the information criteria for a set of models use changeTerms.asrtests with IClikelihood set to REML or full, or use infoCriteria.asreml.

Usage

## S3 method for class 'asrtests'
changeModelOnIC(asrtests.obj, 
                dropFixed = NULL, addFixed = NULL, 
                dropRandom = NULL,  addRandom = NULL, 
                newResidual = NULL, label = "Changed terms", 
                allow.unconverged = TRUE, allow.fixedcorrelation = TRUE,
                checkboundaryonly = FALSE, 
                trace = FALSE, update = TRUE, denDF = "numeric", 
                set.terms = NULL, ignore.suffices = TRUE, 
                bounds = "P", initial.values = NA, 
                which.IC = "AIC", IClikelihood = "REML", 
                fixedDF = NULL, varDF = NULL, 
                bound.exclusions = c("F","B","S","C"), 
          ...)

Arguments

Value

An asrtests.object containing the components (i) asreml.obj, (ii) wald.tab, and (iii) test.summary. The values of the degrees of freedom and the information criteria are differences between those of the changed model and those of the model supplied to changeModelOnIC.

Author(s)

Chris Brien

See Also

as.asrtests, rmboundary.asrtests, testranfix.asrtests, testresidual.asrtests,
newfit.asreml, reparamSigDevn.asrtests, chooseModel.asrtests, changeTerms.asrtests, infoCriteria.asreml

Examples

## Not run: 

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, 
                            label = "Maximal model")
current.asrt <- rmboundary(current.asrt)

# Drop both Row and Column terms
current.asrt <- changeModelOnIC(current.asrt, 
                                dropRandom = "Row + Column", 
                                checkboundaryonly = TRUE,
                                which.IC = "AIC", IClikelihood = "full")
current.asrt <- iterate(current.asrt)
                          
# Add and drop both fixed and random terms
current.asrt <- changeModelOnIC(current.asrt, 
                                addFixed = "vRow", dropFixed = "WithinColPairs", 
                                addRandom = "spl(vRow)", dropRandom = "units", 
                                checkboundaryonly = TRUE,
                                which.IC = "AIC", IClikelihood = "full")
                          
# Replace residual with model without Row autocorrelation
current.asrt <- changeModelOnIC(current.asrt, 
                                newResidual = "Row:ar1(Column)", 
                                label="Row autocorrelation",
                                IClikelihood = "full")


## End(Not run)

Description

The specified terms are simply added or dropped, without testing, from either the fixed or random model and/or the residual (rcov) model replaced. Also, the bounds and/or initial values of some terms can be changed. No hypothesis testing is performed, but a check is made for boundary or singular terms.

A row is added to the test.summary data.frame of the asrtests.object using the supplied label and stating which models have been changed. Information criteria can be included in the row of the test.summary. Convergence in fitting the model is checked and a note included in the action if there was not. All components of the asrtests.object are updated.

To only change the terms based on a comparison of information criteria use changeModelOnIC.asrtests.

Usage

## S3 method for class 'asrtests'
changeTerms(asrtests.obj, 
            dropFixed = NULL, addFixed = NULL, 
            dropRandom = NULL,  addRandom = NULL, 
            newResidual = NULL, label = "Changed terms", 
            allow.unconverged = TRUE, allow.fixedcorrelation = TRUE, 
            checkboundaryonly = FALSE, 
            trace = FALSE, update = TRUE, denDF = "numeric", 
            set.terms = NULL, ignore.suffices = TRUE, 
            bounds = "P", initial.values = NA, 
            IClikelihood = "none", bound.exclusions = c("F","B","S","C"), 
            ...)

Arguments

Value

An asrtests.object containing the components (i) asreml.obj, (ii) wald.tab, and (iii) test.summary.

Author(s)

Chris Brien

References

Kenward, M. G., & Roger, J. H. (1997). Small sample inference for fixed effects from restricted maximum likelihood. Biometrics, 53, 983-997.

See Also

as.asrtests, rmboundary.asrtests, testranfix.asrtests, testresidual.asrtests,
newfit.asreml, reparamSigDevn.asrtests, chooseModel.asrtests,
changeModelOnIC.asrtests, infoCriteria.asreml

Examples

## Not run: 
terms <- "(Date/(Sources * (Type + Species)))"
current.asrt <- changeTerms(current.asrt, addFixed = terms)

current.asrt <- changeTerms(current.asrt, dropFixed = "A + B", denDF = "algebraic")

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(current.asrt)
# Add and drop both fixed and random terms
current.asrt <- changeTerms(current.asrt, 
                            addFixed = "vRow", dropFixed = "WithinColPairs", 
                            addRandom = "spl(vRow)", dropRandom = "units", 
                            checkboundaryonly = TRUE)
# Replace residual with model without Row autocorrelation
current.asrt <- changeTerms(current.asrt, 
                            newResidual = "Row:ar1(Column)", 
                            label="Row autocorrelation")


## End(Not run)

Description

The data collected after imaging had been completed on the 1056 plants in the experiment reported by Atieno et al. (2017). The design employed for the experiment was a split-plot design in which two consecutive carts formed a main plot. The split-plot design assigned 245 genotypes to main plots, the genotypes being unequally replicated 2 or 3 times. Treatments (non-saline, saline) were randomized to the two subplots (carts) within each main plot.

The columns in the data.frame are: Smarthouse, Lane, Position, Zone, Mainplot, Subplot, Replicate, xLane, xPosition, Genotypes, Treatments, Biomass, PlantHeight, SenescenceRank, TotalPods, FilledPods, EmptyPods, SeedNo, TotalSeedWt, SeedWt100.

The columns Smarthouse, Lane and Position uniquely identify the rows of observations. Zones are groups of 4 Lanes, Mainplots are the 44 pairs of consecutive Subplots within each Zone, and a Subplot is a cart containing a single plant. The columns xLane and xPosition are numeric covariates for location within a Smarthouse. Genotypes and Treatments indicate the genotype and treatment that each plant was allocated. The response variables are Biomass, PlantHeight, SenescenceRank, TotalPods, FilledPods, EmptyPods, SeedNo, TotalSeedWt and SeedWt100.

Usage

data(ChickpeaEnd.dat)

Format

A data.frames with 1056 rows by 20 columns.

References

Atieno, J., Li, Y., Langridge, P., Dowling, K., Brien, C., Berger, B., Varshney, R. K., and Sutton, T. (2017). Exploring genetic variation for salinity tolerance in chickpea using image-based phenotyping. Scientific Reports, 7, 1300. doi:10.1038/s41598-017-01211-7

Determines the set of significant terms using p-values and records the tests performed in a data.frame, taking into account the marginality relations of terms.

Description

Using p-values from hypothesis tests, determines the set of significant terms, taking into account the hierarchy or marginality of terms. In particular, a term will not be tested if it is marginal to (or nested in) one that is significant. For example, if A:B is significant, then neither A nor B will be tested. The tests conducted in choosing selected model are listed in a summary data.frame.

Usage

chooseModel(object, ...)

Arguments

Details

chooseModel is the generic function for the chooseModel method. Use methods("chooseModel") to get all the methods for the chooseModel generic.

chooseModel.asrtests is a method for an asrtests.object. It uses testranfix.asrtests to conduct tests to determine the p-values used in the model selection.

chooseModel.data.frame is a method for a data.frame. It uses the p-values stored in the data.frame in the model selection.

Author(s)

Chris Brien

See Also

chooseModel.asrtests, chooseModel.asrtests, changeModelOnIC.asrtests, testranfix.asrtests

Determines and records the set of significant terms using an asrtests.object, taking into account the hierarchy or marginality relations of the terms.

Description

Performs a series of hypothesis tests on a set of fixed and/or random terms taking into account the marginality of terms. In particular, a term will not be tested if it is marginal to (or nested in) one that is significant. For example, if A:B is significant, then neither A nor B will be tested. For a random term, the term is removed from the model fit, any boundary terms are removed using rmboundary.asrtests and a REML likelihood ratio test is performed using REMLRT.asreml. If it is not significant and drop.ran.ns is TRUE, the term is permanently removed from the model. Note that if boundary terms are removed, the reduced model may not be nested in the full model in which case the test is not valid. For fixed terms, the Wald tests are performed and the p-value for the term obtained. If it is not significant and drop.fix.ns is TRUE, the term is permanently removed from the model. A row that records the outcome of a test is added to test.summary for each term that is tested.

Usage

## S3 method for class 'asrtests'
chooseModel(object, terms.marginality=NULL, 
            alpha = 0.05, allow.unconverged = TRUE, 
            allow.fixedcorrelation = TRUE,  
            checkboundaryonly = FALSE, drop.ran.ns=TRUE, 
            positive.zero = FALSE, bound.test.parameters = "none", 
            drop.fix.ns=FALSE, denDF = "numeric",  dDF.na = "none", 
            dDF.values = NULL, trace = FALSE, update = TRUE, 
            set.terms = NULL, ignore.suffices = TRUE, 
            bounds = "P", initial.values = NA, 
            IClikelihood = "none", ...)

Arguments