Version:	4.4.49
Date:	2025-06-10
Title:	Augments 'ASReml-R' in Fitting Mixed Models and Packages Generally in Exploring Prediction Differences
Depends:	R (≥ 3.5.0)
Imports:	dae, devtools, doParallel, dplyr, foreach, ggplot2, graphics, grDevices, methods, nloptr, parallel, qqplotr, RColorBrewer, reshape2, rlang, stats, sticky, stringr, tryCatchLog, utils
Suggests:	emmeans (≥ 1.8.8), lattice, lmerTest, pbkrtest, R.rsp, testthat, tictoc
Enhances:	asreml
VignetteBuilder:	R.rsp
SystemRequirements:	asreml
LazyData:	true
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/.
License:	MIT + file LICENSE
URL:	http://chris.brien.name
BugReports:	https://github.com/briencj/asremlPlus/issues
NeedsCompilation:	no
Packaged:	2025-06-10 08:42:15 UTC; briencj
Author:	Chris Brien [aut, cre]
Maintainer:	Chris Brien <chris.brien@adelaide.edu.au>
Repository:	CRAN
Date/Publication:	2025-06-15 09:10:05 UTC

Augments 'ASReml-R' in Fitting Mixed Models and Packages Generally in Exploring Prediction Differences

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.49

Date: 2025-06-10

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] (ORCID: <https://orcid.org/0000-0003-0581-1817>)

Maintainer: Chris Brien <chris.brien@adelaide.edu.au>

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)

A large data set comprising the end of imaging data from a chick pea experiment conducted in high-throughput greenhouses

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

Description of an LSD frame

Description

A data.frame that stores Least Significant differences (LSDs) for predictions for a fitted model.

Value

A data.frame that can be a component of an alldiffs.object and that contains LSD values and statistics to be used in determining the significance of the pairwise differences. In particular, they are used in calculating halfLeastSignificant limits to be included in a predictions.frame.

Exactly what an LSD.frame contains is determined by the following arguments to functions that return an alldiffs.object: LSDtype, LSDby, LSDstatistic, LSDaccuracy and LSDsupplied. The rownames of the LSD.frame indicate, for each of its rows, for what group of predictions the entries in the row were calculated, this being controlled by the LSDtype and LSDby arguments. The values for all of the LSD arguments are stored as attributes to the alldiffs.object and the predictions and, if present backtransforms, components of the alldiffs.object.

An LSD.frame always has the eight columns c, minimumLSD, meanLSD, maximumLSD, assignedLSD, accuracyLSD, falsePos and falseNeg.

1. c: This gives the number of pairwise comparison of predictions for the combinations of the factor levels given by the row name. If the row name is overall then it is for all predictions.

2. minimumLSD, meanLSD, maximumLSD: These are computed for either overall, factor.combinations, per.prediction or supplied LSD values, as specified by the LSDtype argument. The meanLSD is calculated using the square root of the mean of the variances of set of pairwise differences appropriate to the specific LSDtype argument.

   For overall, the mean, minimum and maximum of the LSDs for all pairwise comparisons are computed.

   If factor.combinations was specified for LSDtype when the LSDs were being calculated, then the LSD.frame contains a row for each combination of the values of the factors and numerics specified by LSDby. The values in a row are calculated from the LSD values for the pairwise differences for each combination of the factors and numerics values, unless there is only one prediction for a combination, when notional LSDs are calculated that are based on the standard error of the prediction multiplied by the square root of two.

   For per.prediction, the minimum, mean and maximum LSD, based, for each prediction, on the LSD values for all pairwise differences involving that prediction are computed.

   For supplied, the LSD.frame is set up based on the setting of LSDby: a single row with name overall if LSDby is NULL or, if LSDby is a vector of factor and numeric names, rows for each observed combinations of the values of the named factors and numerics. The LSDsupplied argument is used to provide the values to be stored in the column assignedLSD.

3. assignedLSD: The assignedLSD column contains the values that are assigned for use in calculating halfLeastSignificant error.intervals. Its contents are determined by LSDstatistic and LSDsupplied arguments. The LSDsupplied argument allows the direct specification of values to be placed in the assignedLSD column of the LSD.frame. The default is to use the values in the meanLSD column.

4. LSDaccuracy: The LSDaccuracy gives an indication of the proportion that the correct LSD for a single predicted.value might deviate from its assignedLSD value. The contents of the accuracyLSD column is controlled by the LSDaccuracy argument.

5. falsePos and falseNeg: These columns 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. Each LSD value in the assignedLSD column is used to determine the significance of pairwise differences that involve predictions for the combination of values given by the row name for the LSD value.

See recalcLSD.alldiffs for more information.

Author(s)

Chris Brien

See Also

recalcLSD.alldiffs, redoErrorIntervals.alldiffs, predictPresent.asreml,
predictPlus.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.diffs <- predictPlus(m1.asr, classify="Nitrogen:Variety", 
                          wald.tab = current.asrt$wald.tab, 
                          tables = "none")
  
## 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)
    #Get predictions
    Var.emm <- emmeans::emmeans(m1.lmer, specs = ~ Nitrogen:Variety)
    Var.preds <- summary(Var.emm)
    ## 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

    #Set up an alldiffs object, which includes overall LSDs
    Var.diffs <- allDifferences(predictions = Var.preds, classify = "Variety:Nitrogen", 
                                     sed = Var.sed, vcov = Var.vcov, tdf = 45)
  }

  if (exists("Var.diffs"))
  {
    ## Use recalcLSD to get LSDs for within Variety differences
    Var.LSD.diffs <- recalcLSD(Var.diffs, 
                               LSDtype = "factor.combinations", LSDby = "Variety")
    print(Var.LSD.diffs$LSD)
  }

Data for an experiment to investigate whether ladybirds transfer aphids

Description

Welham et al. (2015, Example 8.2) describe a three-factor factorial experiment to investigate whether ladybirds transfer fungus to live aphids on plants. The three factors are Host plant (beans, trefoil), infected Cadavers (5, 10, 20), and Ladybird (-, +). A generalized randomized complete-block design is used to assign the three factors to 2 Runs, each of which involves 36 containers with a plant and live aphids. The response to be analyzed is the logit of the proportion of live aphids that were infected.

The columns in the data frame are: ID, Run, Plant, Host, Ladybird, Cadavers, Live, Infected, logitP, Prop. The column ID numbers the observations. Live, Infected, logitP, Prop are response variables.

Usage

data(Ladybird.dat)

Format

A data.frame containing 72 observations of 10 variables.

Author(s)

Chris Brien

Source

Welham, S. J., Gezan, S. A., Clark, S. J., & Mead, A. (2015). Statistical Methods in Biology: Design and Analysis of Experiments and Regression. Boca Raton: Chapman and Hall/CRC..

Data for an experiment to investigate nitrogen response of 3 oats varieties

Description

Yates (1937) describes a split-plot experiment that investigates the effects of three varieties of oats and four levels of Nitrogen fertilizer. The varieties are assigned to the main plots using a randomized complete block design with 6 blocks and the nitrogen levels are randomly assigned to the subplots in each main plot.

The columns in the data frame are: Blocks, Wplots, Subplots, Variety, Nitrogen, xNitrogen, Yield. The column xNitrogen is a numeric version of the factor Nitrogen. The response variable is Yield.

Usage

data(Oats.dat)

Format

A data.frame containing 72 observations of 7 variables.

Author(s)

Chris Brien

Source

Yates, F. (1937). The Design and Analysis of Factorial Experiments. Imperial Bureau of Soil Science, Technical Communication, 35, 1-95.

Calculates the adjusted coefficient of determination for a specified combination of fixed and random terms.

Description

Calculates the adjusted coefficient of determination (R2) that measures the contributions to the total variance exhibited by the observations of a specified combination of fixed and random terms in a fitted linear mixed model.

Note that the adjusted R2 can be negative, which indicates that the contribution of the terms to the total variance is very small relative to the sum of the degrees of freedom of the terms.

Piepho's (2023) method for GLMMs has not been implemented. This function is not available for ASReml-R version 3.

Usage

## S3 method for class 'asreml'
R2adj(asreml.obj, 
      include.which.fixed = ~ ., orthogonalize = "hybrid", 
      include.which.random = NULL, 
      bound.exclusions = c("F","B","S","C"), ...)

Arguments

Details

The method used to compute the adjusted R2 under a linear mixes model (LMM) is that described by Piepho (2023). Here, the method has been extended to allow computation of the adjusted R2 for a subset of the fixed terms. A set of orthogonalized projectors for all of the fixed terms in the model (a set of \mathbf{Q}_i\mathrm{s}) is obtained and the combined contribution of the fixed terms nominated in include.which.fixed is obtained by computing the average semisquared bias, ASSB, for the nominated fixed terms as:

\Sigma_i \{(\mathbf{Q}_i \mathbf{X}\boldsymbol{\beta})^\mathrm{T}\mathbf{Q}_i \mathbf{X}\boldsymbol{\beta} + \textnormal{trace}(\mathbf{X}^\mathrm{T} \mathbf{Q}_i \mathbf{X} \mathrm{var}(\boldsymbol{\beta})) \} / (n - 1)

Of the two methods, eigenmethods is least likely to fail, but it does not establish the marginality between the terms. It is often needed when there is nonorthogonality between terms, such as when there are several linear covariates. It can also be more efficient in these circumstances.

The process can be computationally expensive, particularly for a large data set (500 or more observations) and/or when many terms are to be orthogonalized, particularly if they are not orthogonal.

If the error "Matrix is not idempotent" should occur then, especially if there are many terms, one might try using set.daeTolerance from the dae package to reduce the tolerance used in determining if values are either the same or are zero; it may be necessary to lower the tolerance to as low as 0.001. Also, setting orthogonalize to eigenmethods is worth a try.

In doing the computations, no changes are made to the fitted model, nor is the formula stored in asreml.obj referred to. Instead, the names of the terms referred to are those stored in the coefficients component of the asreml.obj. Use attr(asreml.obj$coefficients$fixed, which = "terms") to access the attribute for fixed terms; substitute random for fixed to see the names of the random terms. For fixed terms. the term names are the same as those in the Wald table produced by wald.asreml, and, for random terms, the same as those in the vparameters component of the asreml.obj. Two asreml formula functions whose terms can differ from their formulation in a model formula are at and str.)

The function estimateV.asreml is used to calculate the variance matrices required in calculating the adjusted R2.

Value

A numeric that is the adjusted R2, expressed as a percentage. It has attributes include.which.fixed, include.which.random and missing.termmatrix (use attr(x, which = "name") to access the attribute name). The missing.termmatrix attribute will be NULL, unless the design matrix could not be obtained for one or more model terms. If is is not NULL, it will be a list of terms whose design matices could not be produced and so are not included in the variance matrix estimate. An NA will be returned for the adjusted R2 if missing.termmatrix is not NULL or a generalized inverse could not be computed for the variance matrix estimate.

Author(s)

Chris Brien

References

Piepho, H.-P. (2023). An adjusted coefficient of determination (R2) for generalized linear mixed models in one go. Biometrical Journal, 65(7), 2200290. doi:10.1002/bimj.202200290.

See Also

asreml, estimateV.asreml.

Examples

## Not run: 
  data(Oats.dat)
  
  current.asr <- asreml(Yield ~ Nitrogen*Variety, 
                        random=~Blocks/Wplots,
                        data=Oats.dat)
  R2.adj.fix <- R2adj.asreml(current.asr)
  R2.adj.ran <- R2adj.asreml(current.asr, 
                             include.which.fixed = NULL, include.which.random = ~ .)
  R2.adj.tot <- R2adj.asreml(current.asr, include.which.random = ~ .)
  R2.adj.tot <- R2adj.asreml(current.asr, include.which.random = ~ Blocks)
  R2.adj.add <- R2adj.asreml(current.asr, include.which.fixed = ~ Nitrogen + Variety)
  R2.adj.int <- R2adj.asreml(current.asr, 
                             include.which.fixed = ~ . - (Nitrogen + Variety))
  R2.adj.int <- R2adj.asreml(current.asr, include.which.fixed = ~ Nitrogen:Variety)

## End(Not run)

Performs a REML ratio test to compare two models.

Description

Extracts the REML log likelihood and the number of variance parameters from two asreml objects. It assumes that the first asreml object corresponds to the null hypothesis and the second asreml object to the alternative hypothesis for the test being conducted. That is, the second asreml object is the result of fitting a model that is a reduced version of the model for the first object. In the case where the reduced model is obtained by setting positively-constrained variance parameters in the full model to zero, the positive.zero argument should be set to TRUE so that the p-value is computed using a mixture of chi-square distributions as described in Self and Liang (1987).

The function checks that the models do not differ in either their fixed or sparse models.

Usage

## S3 method for class 'asreml'
REMLRT(h0.asreml.obj, h1.asreml.obj, 
       positive.zero = FALSE, bound.test.parameters = "none", 
       DF = NULL, bound.exclusions = c("F","B","S","C"), ...)

Arguments

Value

A data.frame containing the log of the likelihood ratio, its degrees of freedom, its p-value and the number of bound parameters in each of the two models being compared.

Note

If DF is not NULL, the supplied value is used. Otherwise DF is determined from the information in h1.asreml.obj and h0.asreml.obj. In this case, the degrees of freedom for the test 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).

The test statistic is calculated as 2(log(REML)_1 - log(REML)_0).

This procedure is only appropriate when the null hypothesis is that (i) all parameters are on the boundary of the parameter space (ii) all parameters are in the interior of the parameter space, or (iii) there are two parameters, one of which is on the boundary and the other is not. Other cases have been discussed by Self and Liang (1987), but are not implemented here.

Author(s)

Chris Brien

References

Self, S.G., and Liang, K-Y. (1987) Asymptotic Properties of Maximum Likelihood Estimators and Likelihood Ratio Tests Under Nonstandard Conditions. Journal of the American Statistical Association, 82, 605-10.

See Also

infoCriteria.asreml, testranfix.asrtests

Examples

## Not run: 
    REMLRT(ICV.max, ICV.red, bound.test.parameters = "onlybound")

## End(Not run)

Data for an experiment to investigate the quality of water runoff over time

Description

This data is from an experiment to investigate the quality of water runoff. However, it has been modified to hide the true identity of the Species and Sources. It is used to provide executable examples of the functions listed under Examples.

Usage

data(WaterRunoff.dat)

Format

A data.frame containing 440 observations of 13 variables.

Author(s)

Chris Brien

Source

Kazemi, F. (pers. comm.)

See Also

chooseModel.asrtests, reparamSigDevn.asrtests,
plotPredictions.data.frame, predictPlus.asreml, predictPresent.asreml

Data for a 1976 experiment to investigate 25 varieties of wheat

Description

The data appears in Gilmour et al. (1995) and is from a field experiment designed to compare the performance of 25 varieties of spring wheat. An analysis of it using asreml is presented by Butler et al. (2023, Section 7.6), although they suggest that it is a barley experiment. It is used in the Wheat vignettes [Enter vignette(package = "asremlPlus")] as an executable example of the use of the asremlPlus to analyse a data set.

The experiment was conducted at Slate Hall Farm, UK, in 1976 and was designed as a balanced lattice square with 6 replicates laid out in a 10 \times 15 rectangular grid. The columns in the data frame are: Rep, Row, Column, WithinColPairs, Variety, yield. The response variable is the grain yield.

Usage

data(Wheat.dat)

Format

A data.frame containing 150 observations of 6 variables.

Author(s)

Chris Brien

Source

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/.

Gilmour, A. R., et al. (1995) Average Information REML: An efficient algorithm for variance parameter estimation in linear mixed models. Biometrics, 51, 1440-1450.

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)
}

Adds, to a supplied model, a spatial model that accounts for local spatial variation.

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)

Uses information criteria to decide whether to add a spatial model to account for local spatial variation.

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)

Adds a row to a test.summary data.frame.

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)

Using supplied predictions and standard errors of pairwise differences or the variance matrix of predictions, forms all pairwise differences between the set of predictions, and p-values for the differences.

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 of an alldiffs object

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)
  }

Applies the angular transformation to proportions.

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)

Applies the modified angular transformation to a vector of counts.

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)
  }

Forms an asrtests object that stores (i) a fitted asreml object, (ii) a pseudo-anova table for the fixed terms and (iii) a history of changes and hypothesis testing used in obtaining the model.

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 that also inherits S3-class list.

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)
  }

Deprecated Functions in the Package asremlPlus

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

The randomly-presented, startup tips.

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 of an asrtests object

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

Uses the parametric bootstrap to calculate the p-value for a REML ratio test to compare two models.

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)

Uses information criteria to decide whether to change an already fitted model.

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, 
                allow.absentDropTerms = FALSE, 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)

Adds and drops terms from one or both of the fixed or random model, replaces the residual (rcov) model with a new model and changes bounds or initial values of terms.

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)

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

Value

A list containing:

1. asrtests.obj: an asrtests.object containing the components (i) asreml.obj, (ii) wald.tab, and (iii) test.summary.;

2. sig.tests: a character vector whose elements are the significant terms amongst those tested.

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

chooseModel, chooseModel.data.frame, as.asrtests, testranfix.asrtests,
testresidual.asrtests, REMLRT.asreml, rmboundary.asrtests, newfit.asreml,
changeModelOnIC.asrtests, changeTerms.asrtests, reparamSigDevn.asrtests

Examples

## Not run: 
data(WaterRunoff.dat)
asreml.options(keep.order = TRUE) #required for asreml-R4 only
current.asr <- asreml(log.Turbidity ~ Benches + (Sources * (Type + Species)) * Date, 
                      random = ~Benches:MainPlots:SubPlots:spl(xDay), 
                      data = WaterRunoff.dat, keep.order = TRUE)
current.asrt <- as.asrtests(current.asr, NULL, NULL)
terms.treat <- c("Sources", "Type", "Species", 
                 "Sources:Type", "Sources:Species")
terms <- sapply(terms.treat, 
                FUN=function(term){paste("Date:",term,sep="")}, 
                simplify=TRUE)
terms <- c("Date", terms)
terms <- unname(terms)
marginality <-  matrix(c(1,0,0,0,0,0, 1,1,0,0,0,0,  1,0,1,0,0,0, 
                         1,0,1,1,0,0, 1,1,1,0,1,0, 1,1,1,1,1,1), nrow=6)
rownames(marginality) <- terms
colnames(marginality) <- terms
choose <- chooseModel(current.asrt, marginality) 
current.asrt <- choose$asrtests.obj
sig.terms <- choose$sig.terms

## End(Not run)

Determines the set of significant terms from results stored in a data.frame, taking into account the marginality relations of terms and recording the tests used in a data.frame.

Description

Uses the p.values from a set of hypothesis tests that are stored in the supplied data.frame to choose a model to describe the effects of the terms corresponding to the p-values, 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 used in choosing the selected model are listed in the data.frame choose.summary.

No change is made to the p.values, the DF and denDF being for information only.

Usage

## S3 method for class 'data.frame'
chooseModel(object, terms=NULL, p.values = "Pr", 
            DF = "Df", denDF = "denDF", omit.DF = FALSE, 
            terms.marginality=NULL, alpha = 0.05, ...)

Arguments

Value

A list containing:

1. choose.summary: a data.frame summarizing the tests carried out in choosing the significant terms; provided omit.DF = FALSE, it has the same columns as a test.summary from an asrtests.object

2. sig.tests: a character vector whose elements are the significant terms amongst those tested.

Author(s)

Chris Brien

See Also

chooseModel, chooseModel.asrtests

Examples

  
  data("Ladybird.dat")
  
  ## Use asreml to get the table of p-values

  ## Not run: 
  m1.asr <- asreml(logitP ~ Host*Cadavers*Ladybird, 
                   random = ~ Run,
                   data = Ladybird.dat)
  current.asrt <- as.asrtests(m1.asr)
  fixed.tab <-  current.asrt$wald.tab
  col.p <- "Pr"
  df = "Df"
  den.df = "denDF"
  
## End(Not run)
  
  ## Use lmeTest to get the table of p-values
  if (requireNamespace("lmerTest", quietly = TRUE) & 
      requireNamespace("emmeans", quietly = TRUE))
  {
    m1.lmer <- lmerTest::lmer(logitP ~ Host*Cadavers*Ladybird + (1|Run),
                              data=Ladybird.dat)
    fixed.tab <- anova(m1.lmer, type = "II")
    col.p <- "Pr(>F)"
    df = "NumDF"
    den.df = "DenDF"
  }
  
  ## Select a model using the table of p-values obtained with either asreml or lmerTest
  if (exists("fixed.tab"))
  {
    term.marg <- dae::marginality(dae::pstructure(~ Host*Cadavers*Ladybird, 
                                                  data = Ladybird.dat))
    chosen <- chooseModel(fixed.tab, p.values = col.p, DF = df, denDF = den.df, 
                          terms.marginality = term.marg)
  }

Uses information criteria to choose the best fitting spatial model for accounting for local spatial variation.

Description

For a response variable measured on a potentially irregular grid of rows and columns of the units, uses information criteria (IC) to decide whether the fit and parsimony of the model fitted to a set of data can be improved by adding, to the fitted model stored in the supplied asrtests.object, one of the following spatial models to account for the local spatial variation: (i) a two-dimensional first-order autocorrelation model, (ii) a two-dimensional tensor-product natural cubic smoothing spline model (TPNCSS), (iii) a two-dimensional tensor-product penalized P-spline model with second-difference penalties (TPPSC2) model, or (iv) a two-dimensional tensor-product penalized linear spline model with first-difference penalties (TPPSL1). The models from which to select can be reduced to a subset of these four models. For each model, a term from the spatial model is only added to the supplied model if the IC of the supplied model is decreased with the addition of that term. If no model improves the IC when a local spatial variation model is added, then the supplied, nonspatial model will be returned. The data can 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 TPPSC2 models, 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.

One or more rows is added to the test.summary data.frame of the asrtests.object, for each section and each spatial model, stating whether or not the new model has been swapped for a model in which the spatial model has been added to 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 any new model.

Usage

## S3 method for class 'asrtests'
chooseSpatialModelOnIC(asrtests.obj, trySpatial = "all", 
                       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), 
                       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", 
                       return.asrts = "best", ...)

Arguments

Details

For each spatial model that is to be fitted, a fitted spatial model is only returned if it improves the fit over and above that achieved with the model fit supplied in the asrtests.obj, because terms in the spatial model are not added unless model fit is improved by their addition as measured by an IC. If return.asrts is all, then this applies to each spatial model specified by trySpatial. To force a spatial model to be fitted use addSpatialModel.asrtests. The model fit supplied in the asrtests.obj should not include terms that will be included in any 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 changing levels in the data frame faster than those of the preceding variables e.g. Row:Column implies that all levels for Column are in consecutive rows of the data.frame that have 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 are used for any 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 chooseSpatialModelOnIC.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 TPPSC2 and TPPSL1 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. Arguments from tpsmmb and changeModelOnIC.asrtests can be supplied in calls to chooseSpatialModelOnIC.asrtests and will be passed on to the relevant function though the ellipses argument (...).

The data for experiment can be divided into sections and an attempt to fit the same spatial model to each is made. The fit may differ for each of the sections, but 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

A list containing four components: (i) asrts, (ii) spatial.IC, (iii) best.spatial.mod, and (iv) best.spatial.IC.

The component asrts itself holds a list of one or more asrtests.objects, either the best overall out of the supplied model and the spatial models, or, for each spatial model, the best out of the supplied model and that spatial model. Each asrtests.object contains the components: (i) asreml.obj, (ii) wald.tab, and (iii) test.summary. If the asrtests.object is the result of fitting a TPPSC2 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.

The spatial.IC component holds a data.frame with summary of the values of the information criteria for the supplied model and those resulting from adding the spatial models to the supplied model. In the case of a spatial correlation model, the information criteria for the selected spatial correlation model is returned. If a spatial model could not be fitted, then all returned values will be NA).

The best.spatial.mod component is a character giving the name of the best spatial model, and best.spatial.AIC gives the value of its AIC.

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,
addSpatialModel.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)

# Choose the best of four models for the local spatial variation
current.asrt <- chooseSpatialModelOnIC(current.asrt, 
                                       row.covar = "cRow", col.covar = "cColumn",
                                       dropRowterm = "Row", dropColterm = "Column",
                                       asreml.option = "grp")

## End(Not run)

Recreates an asreml object so that it is compatible with the currently loaded asreml version.

Description

Recreate an existing asreml object so that it is compatible with the currently loaded asreml version. It involves refitting the model stored in the asreml object.

Usage

## S3 method for class 'asreml'
convAsremlobj(asreml.obj, ...)

Arguments

Value

An asreml object.

Author(s)

Chris Brien

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

newfit.asreml, update.asreml

Examples

## Not run: 
    m1.asr <- convAsremlobj(m1.asr)

## End(Not run)

Converts the effects names for a term stored in the component of an asreml object into a data.frame.

Description

Converts the effects names for a term stored in the component of an asreml object into a data.frame that has a column for each factor and variable in the term. It facilitates adding the effects to the data.frame supplied to asreml for an analysis. This function can only be used with asreml v4.2 or later.

Usage

## S3 method for class 'asreml'
convEffectNames2DataFrame(asreml.obj, term, use = "design.matrix", sep = ":", ...)

Arguments

Value

A data.frame with columns for the factors and variables in term. It includes the attribute effect.names that contains the extracted effects names for the term

Author(s)

Chris Brien

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/.

Examples

## Not run: 
    G.dat <- convEffectNames2DataFrame(m1.asr, term = "Row:Column", use = "G.aom")
    
    G.dat <- lapply(c("at(Smarthouse, 'SW'):Lane:Position", 
                      "at(Smarthouse, 'SE'):Lane:Position"), 
                   function(term, asreml.obj)
                     tmp <- convEffectNames2DataFrame.asreml(asreml.obj, term = term),
                   asreml.obj = m1.asr)
    G.dat <- do.call(rbind, G.dat)

## End(Not run)

Forms the estimated variance, random or residual matrix for the observations from the variance parameter estimates.

Description

Forms the estimated variance (V), random (G) or (R) matrix for the observations, a square symmetric matrix of order equal to the number of observations. The estimates of the variance parameters and the information about the random and residual models for which they were estimated are obtained from the asreml object. This function is not available in ASReml-R version 3.

Usage

## S3 method for class 'asreml'
estimateV(asreml.obj, which.matrix = "V", 
          extra.matrix = NULL, ignore.terms = NULL, fixed.spline.terms = NULL, 
          bound.exclusions = c("F","B","S","C"), ...)

Arguments

Details

The information about the variance parameters in the fitted mixed model are obtained from the G.param and R.param components of the asreml object. The function can deal with the following variance functions in either the random or residual models: id, diag, us, ar1, ar2, ar3, sar,sar2, ma1, ma2, arma, exp, gau, cor, corb and corg. All of these functions, except us, can be combined with either v or h. It will also cope with the following functions in the random model: at, str, spl, dev, grp, fa and rr. Additionally, it can deal with the function dsum in the residual model. For further information see the ASReml-R User Guide Version 4 (Butler et al., 2023).

Value

A matrix containing the estimated variance matrix. It has an attribute missing.termmatrix (use attr(x, which = "missing.termmatrix") to access the atrribute). It will be NULL, unless the design matrix could not be obtained for one or more model terms. If is is not NULL, it will be a list of terms that could not be produced for inclusion in the variance matrix estimate, and NA will be returned for the estimated variance matrix.

Author(s)

Chris Brien

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, simulate.asreml, variofaces.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)
# Form variance matrix based on estimated variance parameters
V <- estimateV(current.asr)


## End(Not run)

Explores the computed LSD values for pairwise differences between predictions.

Description

Given an alldiffs.object with an sed component, the LSDs are calculated for all pairwise comparisons of predictions. It then calculates (i) a table of frequencies of the LSD values, (ii) the distinct values of the LSDs after rounding, (iii) various statistics from the LSD values, (iv) a measure of the accuracy of each of the LSD statistics, (v) the numbers of false positives and false negatives for each of the LSD statistics if pairwise comparisons are based on the LSD statistic, (vi) the accuracy of each statistic in representing the LSD values for each prediction and (vii) a matrix containing the LSD values for comparing each pair of predictions. Histograms of the frequencies can also be produced.

Usage

## S3 method for class 'alldiffs'
exploreLSDs(alldiffs.obj, LSDtype = "overall", LSDby = NULL, 
            LSDaccuracy = "maxAbsDeviation", alpha = 0.05, digits = 3, 
            retain.zeroLSDs = FALSE, 
            zero.tolerance = .Machine$double.eps ^ 0.5, 
            plotHistogram = FALSE, ...)

Arguments

Details

The false positives and negatives are computed by comparing, for each pair of predictions within each levels-combination of the LSDby variables, the significance of the pair difference determined using (i) the true LSD that is computed from the standard error of differences for the pair and (ii) the approximate LSD that is a statistic computed from the true LSDs for all pairwise difference within each levels-combination of the LSDby variables. The number of false positives is the number of pairwise differences for which a difference is declared significant using the approximate LSD, but not using the true LSD. The number of false negatives is the number of pairwise differences for which a difference is declared nonsignificant using the approximate LSD, but significant using the true LSD.

The LSD accuracy for a set of LSDs is a function of the deviations of those LSDs and an LSD statistic calculated from them; the accuracy is expressed as a proportion of the value of the LSD statistic.

Value

A list with components frequencies, distinct.vals, statistics, accuracy, per.pred.accuracy and LSD:

1. frequencies is a data.frame with the frequency distribution of the LSD values;

2. distinct.vals is a list, each component of which contains the distinct values of the LSDs after rounding;

3. statistics is a data.frame with the minimum, quantile10, quantile25, mean, median, quantile75, quantile90, and maximum of the LSD values;

4. accuracy is a data.frame with the accuracies of the minimum, quantile10, quantile25, mean, median, quantile75, quantile90, and maximum of the LSD values with respect to the values from which these statistics are calculated;

5. false.pos is a data.frame with the numbers of false positives for the pairwise comparisons within each levels-combination of the LSDby variables when each of the minimum, quantile10, quantile25, mean, median, quantile75, quantile90, and maximum of the LSD values is used as an approximate LSD in determining the significance of the pairwise differences;

6. false.neg is a data.frame with the numbers of false negatives for the pairwise comparisons within each levels-combination of the LSDby variables when each of the minimum, quantile10, quantile25, mean, median, quantile75, quantile90, and maximum of the LSD values is used as an approximate LSD in determining the significance of the pairwise differences;

7. per.pred.accuracy is a data.frame with the accuracies of the minimum, quantile10, quantile25, mean, median, quantile75, quantile90, and maximum of the LSD values for a set of predictions when these statistics are used to represent the LSDs for the comparisons amongst the set of predictions;

8. LSD is a square matrix containing the LSD values for all pairwise comparisons of the predictions contained in the supplied alldiffs.obj.

In the statistics, accuracy, false.pos and false.neg data.frames, c is the number of pairwise comparisons on which the values in the same row are based. The accuracy measure is specified by the LSDaccuracy argument.

Author(s)

Chris Brien

See Also

asremlPlus-package, plotLSDs.data.frame, plotLSDs.alldiffs,
plotLSDerrors.alldiffs, plotLSDerrors.data.frame, recalcLSD.alldiffs,
pickLSDstatistics.alldiffs, redoErrorIntervals.alldiffs

Examples

data(WaterRunoff.dat)

##Use asreml to get predictions and associated statistics

## Not run: 
asreml.options(keep.order = TRUE) #required for asreml-R4 only
current.asr <- asreml(fixed = pH ~ Benches + (Sources * (Type + Species)), 
                      random = ~ Benches:MainPlots,
                      keep.order=TRUE, data= WaterRunoff.dat)
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(pH ~ Benches + (Sources * (Type + Species)) + 
                              (1|Benches:MainPlots),
                            data=na.omit(WaterRunoff.dat))
  TS.emm <- emmeans::emmeans(m1.lmer, specs = ~ Sources:Type)
  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:Type", 
                             vcov = TS.vcov, tdf = den.df)
  validAlldiffs(TS.diffs)
}  

## Plot p-values for predictions obtained using asreml or lmerTest
if (exists("TS.diffs"))
{
  ##Explore the LSD values for predictions obtained using asreml or lmerTest  
  LSDstat <- exploreLSDs(TS.diffs, LSDtype = "factor.combinations", 
                         LSDby = "Sources")
}

Combines several factors into one in the components of an alldiffs.object

Description

Combines several factors, in the prediction component of object, into one whose levels are the combinations of the used levels of the individual factors. The matching changes are made to the other components and the attributes of the alldiffs.object. If any of the factors to be combined are in LSDby, they are removed from the LSDby, unless the factors to be combined are exactly those in the LSDby. The levels of the factors are combined using fac.combine from the dae package.

Usage

## S3 method for class 'alldiffs'
facCombine(object, factors, order="standard", 
           combine.levels=TRUE, sep="_", level.length = NA, ...)

Arguments

Value

A modified alldiffs.object.

Author(s)

Chris Brien

See Also

as.alldiffs, allDifferences.data.frame, print.alldiffs, sort.alldiffs,
renewClassify.alldiffs; fac.combine in package dae.

Examples

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

  ## Not run: 
  m1.asr <- asreml(logitP ~ Host*Cadavers*Ladybird, 
                   random = ~ Run,
                   data = Ladybird.dat)
  current.asrt <- as.asrtests(m1.asr)
  HCL.pred <- asreml::predict.asreml(m1.asr, classify="Host:Cadavers:Ladybird", 
                                     sed=TRUE)
  HCL.preds <- HCL.pred$pvals
  HCL.sed <- HCL.pred$sed
  HCL.vcov <- NULL
  wald.tab <-  current.asrt$wald.tab
  den.df <- wald.tab[match("Host:Cadavers:Ladybird", rownames(wald.tab)), "denDF"]
  
## 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(logitP ~ Host*Cadavers*Ladybird + (1|Run),
                              data=Ladybird.dat)
    HCL.emm <- emmeans::emmeans(m1.lmer, specs = ~ Host:Cadavers:Ladybird)
    HCL.preds <- summary(HCL.emm)
    den.df <- min(HCL.preds$df)
    ## Modify HCL.preds to be compatible with a predictions.frame
    HCL.preds <- as.predictions.frame(HCL.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
    HCL.vcov <- vcov(HCL.emm)
    HCL.sed <- NULL
  }
  
  ## Use the predictions obtained with either asreml or lmerTest
  if (exists("HCL.preds"))
  {
    ## Form an all.diffs object
    HCL.diffs <- as.alldiffs(predictions = HCL.preds, classify = "Host:Cadavers:Ladybird", 
                             sed = HCL.sed, vcov = HCL.vcov, tdf = den.df)
    
    ## Check the class and validity of the alldiffs object
    is.alldiffs(HCL.diffs)
    validAlldiffs(HCL.diffs)

    ## Combine Cadavers and Ladybird
    HCL.diffs <- facCombine(HCL.diffs, factors = c("Cadavers","Ladybird"))
    
    ## Check the validity of HCL.diffs
    validAlldiffs(HCL.diffs)
  }

Reorders and/or revises the factor levels using the order of old levels in levels.order and the new labels for the levels given in newlabels. The values in levels.order must be unique.

Description

Reorders and revises the levels and labels of a factor, in the prediction component of an alldiffs.object. The values in the levels.order vector should be the same as the levels in the existing factor, but the order can be changed. To revise the levels, specify the new levels in the newlabels vector and these will replace the corresponding value in the levels.order vector. The matching changes are made to the other components and attributes of the alldiffs.object.

Usage

## S3 method for class 'alldiffs'
facRecast(object, factor, levels.order = NULL, newlabels = NULL, ...)

Arguments

Value

A modified alldiffs.object.

Author(s)

Chris Brien

See Also

as.alldiffs, allDifferences.data.frame, print.alldiffs, sort.alldiffs,
facCombine.alldiffs, facRename.alldiffs, renewClassify.alldiffs; fac.recast in package dae.

Examples

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

  ## Not run: 
  m1.asr <- asreml(logitP ~ Host*Cadavers*Ladybird, 
                   random = ~ Run,
                   data = Ladybird.dat)
  current.asrt <- as.asrtests(m1.asr)
  HCL.pred <- asreml::predict.asreml(m1.asr, classify="Host:Cadavers:Ladybird", 
                                     sed=TRUE)
  HCL.preds <- HCL.pred$pvals
  HCL.sed <- HCL.pred$sed
  HCL.vcov <- NULL
  wald.tab <-  current.asrt$wald.tab
  den.df <- wald.tab[match("Host:Cadavers:Ladybird", rownames(wald.tab)), "denDF"]
  
## 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(logitP ~ Host*Cadavers*Ladybird + (1|Run),
                              data=Ladybird.dat)
    HCL.emm <- emmeans::emmeans(m1.lmer, specs = ~ Host:Cadavers:Ladybird)
    HCL.preds <- summary(HCL.emm)
    den.df <- min(HCL.preds$df)
    ## Modify HCL.preds to be compatible with a predictions.frame
    HCL.preds <- as.predictions.frame(HCL.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
    HCL.vcov <- vcov(HCL.emm)
    HCL.sed <- NULL
  }
  
  ## Use the predictions obtained with either asreml or lmerTest
  if (exists("HCL.preds"))
  {
    ## Form an all.diffs object
    HCL.diffs <- allDifferences(predictions = HCL.preds, classify = "Host:Cadavers:Ladybird", 
                                sed = HCL.sed, vcov = HCL.vcov, tdf = den.df)
    
    ## Check the class and validity of the alldiffs object
    is.alldiffs(HCL.diffs)
    validAlldiffs(HCL.diffs)

    ## Recast the Ladybird and Host factors
    HCL.diffs <- facRecast(HCL.diffs, factor = "Ladybird", 
                           newlabels = c("none", "present"))
    HCL.diffs <- facRecast(HCL.diffs, factor = "Ladybird", 
                           levels.order = c("present", "none"), 
                           newlabels = c("yes","no"))
    HCL.diffs <- facRecast.alldiffs(HCL.diffs, factor = "Host", 
                                    levels.order = c("trefoil", "bean"))

    ## Check the validity of HCL.diffs
    validAlldiffs(HCL.diffs)
  }

Renames factors in the prediction component of an alldiffs.object.

Description

Renames factors in the prediction component of an alldiffs.object. These changes are propagated to the other components and attributes of the alldiffs.object.

Usage

## S3 method for class 'alldiffs'
facRename(object, factor.names, newnames,  ...)

Arguments

Value

A modified alldiffs.object.

Author(s)

Chris Brien

See Also

as.alldiffs, allDifferences.data.frame, print.alldiffs, sort.alldiffs,
facCombine.alldiffs, facRecast.alldiffs, renewClassify.alldiffs; fac.recast in package dae.

Examples

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

  ## Not run: 
  m1.asr <- asreml(logitP ~ Host*Cadavers*Ladybird, 
                   random = ~ Run,
                   data = Ladybird.dat)
  current.asrt <- as.asrtests(m1.asr)
  HCL.pred <- asreml::predict.asreml(m1.asr, classify="Host:Cadavers:Ladybird", 
                                     sed=TRUE)
  HCL.preds <- HCL.pred$pvals
  HCL.sed <- HCL.pred$sed
  HCL.vcov <- NULL
  wald.tab <-  current.asrt$wald.tab
  den.df <- wald.tab[match("Host:Cadavers:Ladybird", rownames(wald.tab)), "denDF"]
  
## 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(logitP ~ Host*Cadavers*Ladybird + (1|Run),
                              data=Ladybird.dat)
    HCL.emm <- emmeans::emmeans(m1.lmer, specs = ~ Host:Cadavers:Ladybird)
    HCL.preds <- summary(HCL.emm)
    den.df <- min(HCL.preds$df)
    ## Modify HCL.preds to be compatible with a predictions.frame
    HCL.preds <- as.predictions.frame(HCL.preds, predictions = "emmean", 
                                      se = "SE", interval.type = "CI", 
                                      interval.names = c("lower.CL", "upper.CL"))
    HCL.vcov <- vcov(HCL.emm)
    HCL.sed <- NULL
  }
  
  ## Use the predictions obtained with either asreml or lmerTest
  if (exists("HCL.preds"))
  {
    ## Form an all.diffs object
    HCL.diffs <- allDifferences(predictions = HCL.preds, 
                                classify = "Host:Cadavers:Ladybird", 
                                sed = HCL.sed, vcov = HCL.vcov, tdf = den.df)
    
    ## Check the class and validity of the alldiffs object
    is.alldiffs(HCL.diffs)
    validAlldiffs(HCL.diffs)

    ## Rename Cadavers 
    HCL.diffs <- facRename(HCL.diffs, factor.names = "Cadavers", newnames = "Cadaver.nos")
    
    ## Check the validity of HCL.diffs
    validAlldiffs(HCL.diffs)
  }

Find LSD values that minimize the number of errors in pairwise comparisons of predictions.

Description

Given an alldiffs.object with an sed component, a search is made of a set of equally spaced values between the minimum and maximum values of the LSDs, calculated from the sed component of the alldiffs.object, to identify LSD values that minimize the number of errors made in deciding on the significance of pairs of predicted values stored in the alldiffs.object. If LSDtype is set to overall, a search is made over the range of LSD values for all pairwise comparisons for a single LSD value; if LSDtype is set to factor.combinations, a separate search is made over the LSD values for the set of pairwise comparisons for each factor.combination in order to identify a single value for each set. The number of values used in the search is controlled by the argument nvalues. For each value in the search, the numbers of false positives and false negatives resulting from employing it as the LSD for each set of pairwise comparisons is calculated. A criterion that combines the false positives and negative is calculated using the false.pos.wt, the criterion being the number of false postives times the false.pos.wt plus the number of false negatives. The value chosen for the LSD is the smallest value from amongst those with the minimum value of the criterion and the least number of false positives. A secondary search with 10 equally spaced values is made of the interval below the chosen value and the search value immediately below it to check whether the chosen grid value can be further reduced without changing the value of either its criterion or the number of false positives.

The primary options for changing the numbers of errors associated with the values resulting from the searching is to manipulate the LSDby and/or false.pos.wt arguments.

Usage

## S3 method for class 'alldiffs'
findLSDminerrors(alldiffs.obj, 
                 LSDtype = "overall", LSDby = NULL, 
                 alpha = 0.05, 
                 false.pos.wt = 10, nvalues = 100,
                 retain.zeroLSDs = FALSE, 
                 zero.tolerance = .Machine$double.eps ^ 0.5, 
                 trace = FALSE, ...)

Arguments

Value

A data.frame containing the chosen LSD(s), its(their) numbers of false positives and negatives and the value(s) of the false criterion.

Author(s)

Chris Brien

See Also

asremlPlus-package, exploreLSDs.alldiffs plotLSDs.data.frame, plotLSDs.alldiffs,
plotLSDerrors.alldiffs, plotLSDerrors.data.frame, recalcLSD.alldiffs,
redoErrorIntervals.alldiffs

Examples

data(WaterRunoff.dat)

##Use asreml to get predictions and associated statistics

## Not run: 
asreml.options(keep.order = TRUE) #required for asreml-R4 only
current.asr <- asreml(fixed = pH ~ Benches + (Sources * (Type + Species)), 
                      random = ~ Benches:MainPlots,
                      keep.order=TRUE, data= WaterRunoff.dat)
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(pH ~ Benches + (Sources * (Type + Species)) + 
                              (1|Benches:MainPlots),
                            data=na.omit(WaterRunoff.dat))
  TS.emm <- emmeans::emmeans(m1.lmer, specs = ~ Sources:Type)
  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:Type", 
                             vcov = TS.vcov, tdf = den.df)
  validAlldiffs(TS.diffs)
}  

## Choose LSD values with the minimum mumber of error for pairwise comparisons of 
##   the predictions obtained using asreml or lmerTest
if (exists("TS.diffs"))
{
  ##Pick the LSD values for predictions obtained using asreml or lmerTest  
  minLSD <- findLSDminerrors(TS.diffs)
  TS.diffs <- redoErrorIntervals(TS.diffs, LSDtype = "supplied", LSDsupplied = minLSD["LSD"])
  TS.diffs$LSD
  minLSDs <- findLSDminerrors(TS.diffs, LSDtype = "factor.combinations", 
                              LSDby = "Sources")
  TS.diffs <- redoErrorIntervals(TS.diffs, LSDtype = "supplied", 
                                 LSDby = "Sources", LSDsupplied = minLSDs["LSD"])
  TS.diffs$LSD
}

Finds the version of asreml that is loaded and returns the initial characters in version.

Description

Checks that asreml is loaded and, if it is, returns the first nchar characters of the version that is loaded.

Usage

getASRemlVersionLoaded(nchar = NULL, notloaded.fault = FALSE)

Arguments

Value

A character, being the first nchar characters of the version of asreml that is loaded.

Author(s)

Chris Brien

See Also

loadASRemlVersion.

Examples

## Not run: 
getASRemlVersionLoaded()
## End(Not run)

Gets the formulae from an asreml object.

Description

Gets the formulae nominated in the which argument from the call stored in an asreml object.

Usage

## S3 method for class 'asreml'
getFormulae(asreml.obj, which = c("fixed", "random", "residual"), 
            expanded = FALSE, envir = parent.frame(), ...)

Arguments

Value

A list containing a component with each of the extracted formula(e), the name of a component being the formula that it contains.

Author(s)

Chris Brien

See Also

printFormulae.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)
    getFormulae(current.asr)

## End(Not run)

Gets the entry for a test recorded in the test.summary data.frame of an asrtests.object

Description

Matches the label in the term column of the test.summary data.frame in the supplied asrtests.object and extracts the line for it. It only matches the last occurrence of label.

Usage

## S3 method for class 'asrtests'
getTestEntry(asrtests.obj, label, error.absent = TRUE, ...)

Arguments

Value

A one-line data.frame containing the entry or, error.absent is NULL, NULL.

Author(s)

Chris Brien

See Also

getTestPvalue.asrtests, as.asrtests,
testranfix.asrtests, testswapran.asrtests, testresidual.asrtests,
changeModelOnIC.asrtests, changeTerms.asrtests, chooseModel.asrtests

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)
current.asrt <- rmboundary(current.asrt)
# Test nugget term
current.asrt <- testranfix(current.asrt, "units", positive=TRUE)
getTestEntry(current.asrt, label = "units")

## End(Not run)

Gets the p-value for a test recorded in the test.summary data.frame of an asrtests.object

Description

Matches the label in the term column of the test.summary data.frame in the supplied asrtests.object and extracts its p-value. It only matches the last occurrence of label.

Usage

## S3 method for class 'asrtests'
getTestPvalue(asrtests.obj, label, ...)

Arguments

Value

An numeric containing the p-value. It can be NA, for example when a p-value could not be calculated.

Author(s)

Chris Brien

See Also

getTestEntry.asrtests, as.asrtests,
testranfix.asrtests, testswapran.asrtests, testresidual.asrtests,
changeTerms.asrtests, chooseModel.asrtests

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)
current.asrt <- rmboundary(current.asrt)
# Test nugget term
current.asrt <- testranfix(current.asrt, "units", positive=TRUE)
getTestPvalue(current.asrt, label = "units")

## End(Not run)

Computes AIC and BIC for models.

Description

Computes Akiake and Bayesian (Schwarz) Information Criteria for models. Either the Restricted Maximum likelihood (REML) or the full likelihood (full) can be used. The full likelihood, evaluated using REML estimates is used when it is desired to compare models that differ in their fixed models.

Usage

## S3 method for class 'asreml'
infoCriteria(object, DF = NULL, 
            bound.exclusions = c("F","B","S","C"), 
            IClikelihood = "REML", fixedDF = NULL, varDF = NULL, ...)
## S3 method for class 'list'
infoCriteria(object, bound.exclusions = c("F","B","S","C"), 
            IClikelihood = "REML", fixedDF = NULL, varDF = NULL, ...)

Arguments