Multiple Comparison Criteria for ASReml-R Predicted Values
compare.RdComputes a pairwise comparison criterion (HSD, LSD, or Bonferroni-corrected
LSD) for predicted values obtained from an ASReml-R model, optionally
within groups defined by a subset of the factors in term.
For each group the function computes the average SED across all \(\binom{n_g}{2}\) pairs of predictions:
$$ \overline{\text{SED}} = \frac{1}{\binom{n_g}{2}} \sum_{i < j} \sqrt{V_{ii} + V_{jj} - 2V_{ij}} $$
where \(V_{ij}\) is the \((i,j)\) element of the prediction error variance-covariance matrix. The comparison criterion is then:
"HSD"(Tukey)\(\displaystyle\frac{\overline{\text{SED}}}{\sqrt{2}} \times q_{1-\alpha}(n_g,\, \nu)\), where \(q\) is the Studentised range quantile with \(n_g\) means and \(\nu = \texttt{model\$nedf}\) error df. Two means are declared significantly different if their absolute difference exceeds the HSD.
"LSD"(unadjusted)\(\overline{\text{SED}} \times t_{\alpha/2,\,\nu}\). No multiplicity correction; liberal for large numbers of comparisons.
"Bonferroni"\(\overline{\text{SED}} \times t_{\alpha/(2m),\,\nu}\), where \(m = \binom{n_g}{2}\) is the number of pairs. Controls the family-wise error rate.
Usage
compare(
model,
term,
by = NULL,
type = c("HSD", "LSD", "Bonferroni"),
pev = TRUE,
alpha = 0.05,
...
)Arguments
- model
An ASReml-R V4 model object.
- term
Character string specifying the
classifyterm, e.g."Treatment:Genotype". Predictions are obtained viaasreml::predict.asreml().- by
Grouping specification. Comparisons are performed independently within each level of the group. Accepted forms:
NULL(default)All predictions form one group.
- Single string
"Site" Split by the
Sitefactor interm.- Colon string
"Site:Treatment" Combine
SiteandTreatmentinto a single grouping variable and split by their interaction.- Character vector
c("Site", "Treatment") Equivalent to the colon form above.
All
byvariables must be present interm. At least one variable intermmust remain for comparisons to be made.- type
Comparison type. One of
"HSD"(Tukey, default),"LSD"(unadjusted), or"Bonferroni".- pev
Logical. If
TRUE(default) the prediction error variance (PEV) frompredict()$vcovis used. IfFALSEandtermis in the random part of the model, the posterior variance \(\bm{G} - \text{PEV}\) is used instead.- alpha
Significance level. Default
0.05.- ...
Additional arguments forwarded to
asreml::predict.asreml().
Value
A data frame with one row per non-missing predicted value, containing:
- Factor columns
All factor variables from
term.predicted.valueThe BLUE or BLUP.
std.errorPrediction standard error.
{type}The comparison criterion (constant within each group). Two predicted values in the same group are significantly different if their absolute difference exceeds this value.
avsedThe average SED for the group (constant within group).
See also
waldTest() for pairwise p-values, asreml::predict.asreml()
Examples
if (FALSE) { # \dontrun{
## HSD for Treatment:Genotype, one group per Site
res <- compare(model,
term = "Treatment:Site:Genotype",
by = "Site",
type = "HSD")
## LSD computed within each Treatment x Site combination
res2 <- compare(model,
term = "Treatment:Site:Genotype",
by = c("Treatment", "Site"),
type = "LSD")
## Equivalently using colon notation
res3 <- compare(model,
term = "Treatment:Site:Genotype",
by = "Treatment:Site",
type = "LSD")
## Bonferroni-corrected LSD, single global group
res4 <- compare(model,
term = "Treatment:Genotype",
type = "Bonferroni",
alpha = 0.05)
## Posterior variance (BLUPs)
res5 <- compare(model,
term = "Treatment:Genotype",
by = "Treatment",
pev = FALSE)
} # }