Nothing
R/ammistability.R
In ammistability: Additive Main Effects and Multiplicative Interaction Model Stability Parameters
Defines functions
ammistability
Documented in
ammistability
#' Estimate multiple AMMI model Stability Parameters
#'
#' \code{ammistability} computes multiple stability parameters from an AMMI
#' model. Further, the corresponding Simultaneous Selection Indices for Yield
#' and Stability (SSI) are also calculated according to the argument
#' \code{ssi.method}. From the results, correlation between the computed indices
#' will also be computed. The resulting correlation matrices will be plotted as
#' correlograms. For visual comparisons of ranks of genotypes for different
#' indices, slopegraphs and heatmaps will also be generated by this function.
#'
#' \code{ammistability} computes the following stability parameters from an AMMI
#' model.
#'
#' \describe{ \item{\link[ammistability:AMGE.AMMI]{Sum Across Environments of
#' GEI Modelled by AMMI
#' (AMGE)}}{\insertCite{sneller_repeatability_1997;textual}{ammistability}}
#' \item{\link[ammistability:ASI.AMMI]{AMMI Stability Index
#' (ASI)}}{\insertCite{jambhulkar_ammi_2014,jambhulkar_genotype_2015,jambhulkar_stability_2017;textual}{ammistability}}
#' \item{\link[ammistability:MASV.AMMI]{AMMI Stability Value
#' (ASV)}}{\insertCite{purchase_parametric_1997,purchase_use_1999,purchase_genotype_2000;textual}{ammistability}}
#' \item{\link[ammistability:ASTAB.AMMI]{AMMI Based Stability Parameter
#' (ASTAB)}}{\insertCite{rao_use_2005;textual}{ammistability}}
#' \item{\link[ammistability:AVAMGE.AMMI]{Sum Across Environments of Absolute
#' Value of GEI Modelled by AMMI
#' (AVAMGE)}}{\insertCite{zali_evaluation_2012;textual}{ammistability}}
#' \item{\link[ammistability:DA.AMMI]{Annicchiarico's D Parameter
#' (DA)}}{\insertCite{annicchiarico_joint_1997;textual}{ammistability}}
#' \item{\link[ammistability:DZ.AMMI]{Zhang's D Parameter
#' (DZ)}}{\insertCite{zhang_analysis_1998;textual}{ammistability}}
#' \item{\link[ammistability:EV.AMMI]{Averages of the Squared Eigenvector Values
#' (EV)}}{\insertCite{zobel_stress_1994;textual}{ammistability}}
#' \item{\link[ammistability:FA.AMMI]{Stability Measure Based on Fitted AMMI
#' Model (FA)}}{\insertCite{raju_study_2002;textual}{ammistability}}
#' \item{\link[ammistability:MASI.AMMI]{Modified AMMI Stability Index
#' (MASI)}}{\insertCite{ajay_modified_2018;textual}{ammistability}}
#' \item{\link[ammistability:MASV.AMMI]{Modified AMMI Stability Value
#' (MASV)}}{\insertCite{zali_evaluation_2012,ajay2019ammistability;textual}{ammistability}}
#' \item{\link[ammistability:SIPC.AMMI]{Sums of the Absolute Value of the IPC
#' Scores
#' (SIPC)}}{\insertCite{sneller_repeatability_1997;textual}{ammistability}}
#' \item{\link[ammistability:ZA.AMMI]{Absolute Value of the Relative
#' Contribution of IPCs to the Interaction
#' (Za)}}{\insertCite{zali_evaluation_2012;textual}{ammistability}} }
#'
#' @inheritParams MASV.AMMI
#' @param AMGE If \code{TRUE}, computes AMGE (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param ASI If \code{TRUE}, computes ASI (see \strong{Details}). \code{n = 2}
#' will be used in this case. Default is \code{TRUE}.
#' @param ASV If \code{TRUE}, computes ASV (see \strong{Details}). \code{n = 2}
#' will be used in this case. Default is \code{TRUE}.
#' @param ASTAB If \code{TRUE}, computes ASTAB (see \strong{Details}). Default
#' is \code{TRUE}.
#' @param AVAMGE If \code{TRUE}, computes AVAMGE (see \strong{Details}). Default
#' is \code{TRUE}.
#' @param DA If \code{TRUE}, computes DA (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param DZ If \code{TRUE}, computes DZ (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param EV If \code{TRUE}, computes EV (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param FA If \code{TRUE}, computes FA (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param MASI If \code{TRUE}, computes MASI (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param MASV If \code{TRUE}, computes MASV (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param SIPC If \code{TRUE}, computes SIPC (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param ZA If \code{TRUE}, computes ZA (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param force.grouping If \code{TRUE}, genotypes will be considered as a
#' grouping variable for plotting the slopegraphs. (Each genotype will be
#' represented by a different colour in the slopegraphs). Default is
#' \code{TRUE}.
#' @param line.size Size of lines plotted in the slopegraphs. Must be numeric.
#' @param line.alpha Transparency of lines plotted in the slopegraphs. Must be
#' numeric.
#' @param line.col Default is \code{TRUE}. Overrides colouring by
#' \code{force.grouping} argument.
#' @param point.size Size of points plotted in the slopegraphs. Must be numeric.
#' @param point.alpha Transparency of points plotted in the slopegraphs. Must be
#' numeric.
#' @param point.col Default is \code{TRUE}. Overrides colouring by
#' \code{force.grouping} argument.
#' @param text.size Size of text annotations plotted in the slopegraphs. Must be
#' numeric.
#'
#' @return A list with the following components: \item{Details}{A data frame
#' indicating the stability parameters computed and the method used for
#' computing the SSI.} \item{Stability Parameters}{A data frame of computed
#' stability parameters.} \item{Simultaneous Selection Indices}{A data frame
#' of computed SSIs.} \item{SP Correlation}{A data frame of correlation
#' between stability parameters.} \item{SSI Correlation}{A data frame of
#' correlation between SSIs.} \item{SP and SSI Correlation}{A data frame of
#' correlation between stability parameters and SSIs.} \item{SP
#' Correlogram}{Correlogram of stability parameters.} \item{SSI
#' Correlogram}{Correlogram of SSIs.} \item{SP and SSI
#' Correlogram}{Correlogram of stability parameters and SSIs.} \item{SP
#' Slopegraph}{Slopegraph of stability parameter ranks.} \item{SSI
#' Slopegraph}{Slopegraph of SSI ranks.} \item{SP Heatmap}{Heatmap of
#' stability parameter ranks.} \item{SSI Heatmap}{Heatmap of SSI ranks.}
#'
#' @importFrom methods is
#' @importFrom stats cor
#' @importFrom stats reorder
#' @import ggcorrplot
#' @importFrom reshape2 melt
#' @export
#'
#' @seealso \code{\link[agricolae]{AMMI}},
#' \code{\link[ammistability]{AMGE.AMMI}},
#' \code{\link[ammistability]{ASI.AMMI}},
#' \code{\link[ammistability]{ASTAB.AMMI}},
#' \code{\link[ammistability]{AMGE.AMMI}},
#' \code{\link[ammistability]{DA.AMMI}}, \code{\link[ammistability]{DZ.AMMI}},
#' \code{\link[ammistability]{EV.AMMI}}, \code{\link[ammistability]{FA.AMMI}},
#' \code{\link[ammistability]{MASV.AMMI}},
#' \code{\link[ammistability]{SIPC.AMMI}},
#' \code{\link[ammistability]{ZA.AMMI}}, \code{\link[ammistability]{SSI}}
#'
#' @references
#'
#' \insertAllCited{}
#'
#' @examples
#' library(agricolae)
#' data(plrv)
#'
#' # AMMI model
#' model <- with(plrv, AMMI(Locality, Genotype, Rep, Yield, console = FALSE))
#'
#' ammistability(model, AMGE = TRUE, ASI = FALSE, ASV = TRUE, ASTAB = FALSE,
#' AVAMGE = FALSE, DA = FALSE, DZ = FALSE, EV = TRUE,
#' FA = FALSE, MASI = FALSE, MASV = TRUE, SIPC = TRUE,
#' ZA = FALSE)
ammistability <- function(model, n, alpha = 0.05,
ssi.method = c("farshadfar", "rao"), a = 1,
AMGE = TRUE, ASI = TRUE, ASV = TRUE, ASTAB = TRUE,
AVAMGE = TRUE, DA = TRUE, DZ = TRUE, EV = TRUE,
FA = TRUE, MASI = TRUE, MASV = TRUE, SIPC = TRUE,
ZA = TRUE, force.grouping = TRUE,
line.size = 1, line.alpha = 0.5, line.col = NULL,
point.size = 1, point.alpha = 0.5, point.col = NULL,
text.size = 2) {
# Check model class
if (!is(model, "AMMI")) {
stop('"model" is not of class "AMMI"')
}
# Check alpha value
if (!(0 < alpha && alpha < 1)) {
stop('"alpha" should be between 0 and 1 (0 < alpha < 1)')
}
# Find number of significant IPCs according to F test
if (missing(n) || is.null(n)) {
n <- sum(model$analysis$Pr.F <= alpha, na.rm = TRUE)
}
# Check for n
if (n %% 1 != 0 && length(n) != 1) {
stop('"n" is not an integer vector of unit length')
}
ssi.method <- match.arg(ssi.method)
sps <- data.frame(sp = c("AMGE", "ASI", "ASV", "ASTAB", "AVAMGE", "DA", "DZ",
"EV", "FA", "MASI", "MASV", "SIPC", "ZA"),
logical = c(AMGE, ASI, ASV, ASTAB, AVAMGE, DA, DZ,
EV, FA, MASI, MASV, SIPC, ZA))
if (nrow(sps[sps$logical, ]) < 2) {
stop('At least two stability parameters should be selected for computation')
}
means <- ASI.AMMI(model)
means$genotype <- rownames(means)
rownames(means) <- NULL
means <- means[, c("genotype", "means")]
SP <- means
SSI <- means
if (AMGE == TRUE) {
res <- AMGE.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "AMGE")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "AMGE"
rm(res)
}
if (ASI == TRUE) {
res <- MASI.AMMI(model, n = 2, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
colnames(res)[colnames(res) == "MASI"] <- "ASI"
SP <- merge(SP, res[, c("genotype", "ASI")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "ASI"
rm(res)
}
if (ASV == TRUE) {
res <- MASV.AMMI(model, n = 2, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
colnames(res)[colnames(res) == "MASV"] <- "ASV"
SP <- merge(SP, res[, c("genotype", "ASV")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "ASV"
rm(res)
}
if (ASTAB == TRUE) {
res <- ASTAB.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "ASTAB")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "ASTAB"
rm(res)
}
if (AVAMGE == TRUE) {
res <- AVAMGE.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "AVAMGE")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "AVAMGE"
rm(res)
}
if (DA == TRUE) {
res <- DA.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "DA")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "DA"
rm(res)
}
if (DZ == TRUE) {
res <- DZ.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "DZ")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "DZ"
rm(res)
}
if (EV == TRUE) {
res <- EV.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "EV")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "EV"
rm(res)
}
if (FA == TRUE) {
res <- FA.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "FA")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "FA"
rm(res)
}
if (MASI == TRUE) {
res <- MASI.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "MASI")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "MASI"
rm(res)
}
if (MASV == TRUE) {
res <- MASV.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "MASV")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "MASV"
rm(res)
}
if (SIPC == TRUE) {
res <- SIPC.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "SIPC")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "SIPC"
rm(res)
}
if (ZA == TRUE) {
res <- ZA.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "Za")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "Za"
rm(res)
}
# Correlations
#-----------------------------------------------------------------------------
# SP corr
del <- names(SP) %in% c("genotype", "means")
spcorr <- cor(SP[, !del])
sppmat <- ggcorrplot::cor_pmat(SP[, !del])
spcorrdf <- formatC(round(spcorr, 2), digits = 2, format = "f")
spcorrdf[] <- paste0(spcorrdf,
ifelse(sppmat < .01, "**", ifelse(sppmat < .05, "*", "")),
sep = "")
spcorrdf[upper.tri(spcorrdf)] <- NA
spcorrdf <- data.frame(spcorrdf, stringsAsFactors = FALSE)
spcorrg <- ggcorrplot(spcorr, hc.order = FALSE, type = "lower",
outline.color = "white", p.mat = sppmat,
ggtheme = theme_bw, show.diag = TRUE, lab = TRUE)
spcorrg <- spcorrg +
ggtitle("Correlation between different\nAMMI stability parameters")
spcorrg <- spcorrg +
theme(axis.text = element_text(colour = "black"))
# SSI corr
del <- names(SSI) %in% c("genotype", "means")
ssicorr <- cor(SSI[, !del])
ssipmat <- ggcorrplot::cor_pmat(SSI[, !del])
ssicorrdf <- formatC(round(ssicorr, 2), digits = 2, format = "f")
ssicorrdf[] <- paste0(ssicorrdf,
ifelse(ssipmat < .01, "**",
ifelse(ssipmat < .05, "*", "")),
sep = "")
ssicorrdf[upper.tri(ssicorrdf)] <- NA
ssicorrdf <- data.frame(ssicorrdf, stringsAsFactors = FALSE)
ssicorrg <- ggcorrplot(ssicorr, hc.order = FALSE, type = "lower",
outline.color = "white", p.mat = ssipmat,
ggtheme = theme_bw, show.diag = TRUE, lab = TRUE)
ssicorrg <- ssicorrg +
ggtitle("Correlation between simultaneous selection indices\nfrom different AMMI stability parameters")
ssicorrg <- ssicorrg +
theme(axis.text = element_text(colour = "black"))
# SP SSI corr
colnames(SSI) <- paste(colnames(SSI), "_SSI", sep = "")
colnames(SSI)[colnames(SSI) == "genotype_SSI"] <- "genotype"
colnames(SSI)[colnames(SSI) == "means_SSI"] <- "means"
spssi <- merge(SP, SSI[, names(SSI)[!(names(SSI) %in% "means")]],
by = "genotype")
del <- names(spssi) %in% c("genotype", "means")
spssicorr <- cor(spssi[, !del])
spssipmat <- ggcorrplot::cor_pmat(spssi[, !del])
spssicorrdf <- formatC(round(spssicorr, 2), digits = 2, format = "f")
spssicorrdf[] <- paste0(spssicorrdf,
ifelse(spssipmat < .01, "**",
ifelse(spssipmat < .05, "*", "")),
sep = "")
spssicorrdf[upper.tri(spssicorrdf)] <- NA
spssicorrdf <- data.frame(spssicorrdf, stringsAsFactors = FALSE)
spssicorrg <- ggcorrplot(spssicorr, hc.order = FALSE, type = "lower",
outline.color = "white", p.mat = spssipmat,
ggtheme = theme_bw, show.diag = TRUE, lab = TRUE,
title = )
spssicorrg <- spssicorrg +
ggtitle("Correlation between different AMMI stability parameters\nand corresponding simultaneous selection indices")
spssicorrg <- spssicorrg +
theme(axis.text = element_text(colour = "black"))
# Bump plots/ slopegraph
#-----------------------------------------------------------------------------
SPrank <- rankdf(df = SP, decreasing = "means",
increasing = names(SP)[!(names(SP) %in% c("genotype",
"means"))])
SPslopeg <- rankslopegraph(df = SPrank, names = "genotype",
force.grouping = TRUE, point.size = point.size,
point.alpha = point.alpha, line.size = line.size,
line.alpha = line.alpha, text.size = text.size,
legend.position = "none") +
ggtitle("Slopegraph of ranks of mean yields and AMMI stability parameters")
SPslopeg <- SPslopeg +
theme(axis.text = element_text(colour = "black"))
SSIrank <- rankdf(df = SSI, decreasing = "means",
increasing = names(SSI)[!(names(SSI) %in% c("genotype",
"means"))])
SSIslopeg <- rankslopegraph(df = SSIrank, names = "genotype",
force.grouping = TRUE, point.size = point.size,
point.alpha = point.alpha, line.size = line.size,
line.alpha = line.alpha, text.size = text.size,
legend.position = "none") +
ggtitle("Slopegraph of ranks of mean yields and simultaneous selction indices")
SSIslopeg <- SSIslopeg +
theme(axis.text = element_text(colour = "black"))
# heatmap
#-----------------------------------------------------------------------------
SPrank$genotype <- reorder(SPrank$genotype, SPrank$means)
SPrankmelt <- melt(SPrank[order(SPrank$means),], id.vars = "genotype")
SPheat <- ggplot(SPrankmelt, aes(x = variable, y = genotype)) +
geom_tile(aes(fill = value)) +
scale_x_discrete(expand = c(0, 0)) +
scale_y_discrete(expand = c(0, 0)) +
scale_fill_gradient(low = "black", high = "green") +
ylab("Genotype") + labs(fill = "Rank") +
theme_bw() +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
axis.text = element_text(colour = "black"),
axis.title.x = element_blank(),
legend.position = "bottom")
SSIrank$genotype <- reorder(SSIrank$genotype, SSIrank$means)
SSIrankmelt <- melt(SSIrank[order(SSIrank$means),], id.vars = "genotype")
SSIheat <- ggplot(SSIrankmelt, aes(x = variable, y = genotype)) +
geom_tile(aes(fill = value)) +
scale_x_discrete(expand = c(0, 0)) +
scale_y_discrete(expand = c(0, 0)) +
scale_fill_gradient(low = "black", high = "green") +
ylab("Genotype") + labs(fill = "Rank") +
theme_bw() +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
axis.text = element_text(colour = "black"),
axis.title.x = element_blank(),
legend.position = "bottom")
details <- list(`Stability parameters estimated` = as.character(sps[sps$logical, ]$sp),
`SSI method` = ifelse(ssi.method == "rao",
paste("Rao and Prabhakaran (2005)",
"with a =", a),
"Farshadfar (2008)"))
out <- list(Details = details,
`Stability Parameters` = SP,
`Simultaneous Selection Indices` = SSI,
`SP Correlation` = spcorrdf,
`SSI Correlation` = ssicorrdf,
`SP and SSI Correlation` = spssicorrdf,
`SP Correlogram` = spcorrg,
`SSI Correlogram` = ssicorrg,
`SP and SSI Correlogram` = spssicorrg,
`SP Slopegraph` = SPslopeg,
`SSI Slopegraph` = SSIslopeg,
`SP Heatmap` = SPheat,
`SSI Heatmap` = SSIheat)
return(out)
}
Try the ammistability package in your browser
Any scripts or data that you put into this service are public.
ammistability documentation built on May 31, 2023, 7:44 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions ammistability
Documented in ammistability
#' Estimate multiple AMMI model Stability Parameters
#'
#' \code{ammistability} computes multiple stability parameters from an AMMI
#' model. Further, the corresponding Simultaneous Selection Indices for Yield
#' and Stability (SSI) are also calculated according to the argument
#' \code{ssi.method}. From the results, correlation between the computed indices
#' will also be computed. The resulting correlation matrices will be plotted as
#' correlograms. For visual comparisons of ranks of genotypes for different
#' indices, slopegraphs and heatmaps will also be generated by this function.
#'
#' \code{ammistability} computes the following stability parameters from an AMMI
#' model.
#'
#' \describe{ \item{\link[ammistability:AMGE.AMMI]{Sum Across Environments of
#' GEI Modelled by AMMI
#' (AMGE)}}{\insertCite{sneller_repeatability_1997;textual}{ammistability}}
#' \item{\link[ammistability:ASI.AMMI]{AMMI Stability Index
#' (ASI)}}{\insertCite{jambhulkar_ammi_2014,jambhulkar_genotype_2015,jambhulkar_stability_2017;textual}{ammistability}}
#' \item{\link[ammistability:MASV.AMMI]{AMMI Stability Value
#' (ASV)}}{\insertCite{purchase_parametric_1997,purchase_use_1999,purchase_genotype_2000;textual}{ammistability}}
#' \item{\link[ammistability:ASTAB.AMMI]{AMMI Based Stability Parameter
#' (ASTAB)}}{\insertCite{rao_use_2005;textual}{ammistability}}
#' \item{\link[ammistability:AVAMGE.AMMI]{Sum Across Environments of Absolute
#' Value of GEI Modelled by AMMI
#' (AVAMGE)}}{\insertCite{zali_evaluation_2012;textual}{ammistability}}
#' \item{\link[ammistability:DA.AMMI]{Annicchiarico's D Parameter
#' (DA)}}{\insertCite{annicchiarico_joint_1997;textual}{ammistability}}
#' \item{\link[ammistability:DZ.AMMI]{Zhang's D Parameter
#' (DZ)}}{\insertCite{zhang_analysis_1998;textual}{ammistability}}
#' \item{\link[ammistability:EV.AMMI]{Averages of the Squared Eigenvector Values
#' (EV)}}{\insertCite{zobel_stress_1994;textual}{ammistability}}
#' \item{\link[ammistability:FA.AMMI]{Stability Measure Based on Fitted AMMI
#' Model (FA)}}{\insertCite{raju_study_2002;textual}{ammistability}}
#' \item{\link[ammistability:MASI.AMMI]{Modified AMMI Stability Index
#' (MASI)}}{\insertCite{ajay_modified_2018;textual}{ammistability}}
#' \item{\link[ammistability:MASV.AMMI]{Modified AMMI Stability Value
#' (MASV)}}{\insertCite{zali_evaluation_2012,ajay2019ammistability;textual}{ammistability}}
#' \item{\link[ammistability:SIPC.AMMI]{Sums of the Absolute Value of the IPC
#' Scores
#' (SIPC)}}{\insertCite{sneller_repeatability_1997;textual}{ammistability}}
#' \item{\link[ammistability:ZA.AMMI]{Absolute Value of the Relative
#' Contribution of IPCs to the Interaction
#' (Za)}}{\insertCite{zali_evaluation_2012;textual}{ammistability}} }
#'
#' @inheritParams MASV.AMMI
#' @param AMGE If \code{TRUE}, computes AMGE (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param ASI If \code{TRUE}, computes ASI (see \strong{Details}). \code{n = 2}
#' will be used in this case. Default is \code{TRUE}.
#' @param ASV If \code{TRUE}, computes ASV (see \strong{Details}). \code{n = 2}
#' will be used in this case. Default is \code{TRUE}.
#' @param ASTAB If \code{TRUE}, computes ASTAB (see \strong{Details}). Default
#' is \code{TRUE}.
#' @param AVAMGE If \code{TRUE}, computes AVAMGE (see \strong{Details}). Default
#' is \code{TRUE}.
#' @param DA If \code{TRUE}, computes DA (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param DZ If \code{TRUE}, computes DZ (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param EV If \code{TRUE}, computes EV (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param FA If \code{TRUE}, computes FA (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param MASI If \code{TRUE}, computes MASI (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param MASV If \code{TRUE}, computes MASV (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param SIPC If \code{TRUE}, computes SIPC (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param ZA If \code{TRUE}, computes ZA (see \strong{Details}). Default is
#' \code{TRUE}.
#' @param force.grouping If \code{TRUE}, genotypes will be considered as a
#' grouping variable for plotting the slopegraphs. (Each genotype will be
#' represented by a different colour in the slopegraphs). Default is
#' \code{TRUE}.
#' @param line.size Size of lines plotted in the slopegraphs. Must be numeric.
#' @param line.alpha Transparency of lines plotted in the slopegraphs. Must be
#' numeric.
#' @param line.col Default is \code{TRUE}. Overrides colouring by
#' \code{force.grouping} argument.
#' @param point.size Size of points plotted in the slopegraphs. Must be numeric.
#' @param point.alpha Transparency of points plotted in the slopegraphs. Must be
#' numeric.
#' @param point.col Default is \code{TRUE}. Overrides colouring by
#' \code{force.grouping} argument.
#' @param text.size Size of text annotations plotted in the slopegraphs. Must be
#' numeric.
#'
#' @return A list with the following components: \item{Details}{A data frame
#' indicating the stability parameters computed and the method used for
#' computing the SSI.} \item{Stability Parameters}{A data frame of computed
#' stability parameters.} \item{Simultaneous Selection Indices}{A data frame
#' of computed SSIs.} \item{SP Correlation}{A data frame of correlation
#' between stability parameters.} \item{SSI Correlation}{A data frame of
#' correlation between SSIs.} \item{SP and SSI Correlation}{A data frame of
#' correlation between stability parameters and SSIs.} \item{SP
#' Correlogram}{Correlogram of stability parameters.} \item{SSI
#' Correlogram}{Correlogram of SSIs.} \item{SP and SSI
#' Correlogram}{Correlogram of stability parameters and SSIs.} \item{SP
#' Slopegraph}{Slopegraph of stability parameter ranks.} \item{SSI
#' Slopegraph}{Slopegraph of SSI ranks.} \item{SP Heatmap}{Heatmap of
#' stability parameter ranks.} \item{SSI Heatmap}{Heatmap of SSI ranks.}
#'
#' @importFrom methods is
#' @importFrom stats cor
#' @importFrom stats reorder
#' @import ggcorrplot
#' @importFrom reshape2 melt
#' @export
#'
#' @seealso \code{\link[agricolae]{AMMI}},
#' \code{\link[ammistability]{AMGE.AMMI}},
#' \code{\link[ammistability]{ASI.AMMI}},
#' \code{\link[ammistability]{ASTAB.AMMI}},
#' \code{\link[ammistability]{AMGE.AMMI}},
#' \code{\link[ammistability]{DA.AMMI}}, \code{\link[ammistability]{DZ.AMMI}},
#' \code{\link[ammistability]{EV.AMMI}}, \code{\link[ammistability]{FA.AMMI}},
#' \code{\link[ammistability]{MASV.AMMI}},
#' \code{\link[ammistability]{SIPC.AMMI}},
#' \code{\link[ammistability]{ZA.AMMI}}, \code{\link[ammistability]{SSI}}
#'
#' @references
#'
#' \insertAllCited{}
#'
#' @examples
#' library(agricolae)
#' data(plrv)
#'
#' # AMMI model
#' model <- with(plrv, AMMI(Locality, Genotype, Rep, Yield, console = FALSE))
#'
#' ammistability(model, AMGE = TRUE, ASI = FALSE, ASV = TRUE, ASTAB = FALSE,
#' AVAMGE = FALSE, DA = FALSE, DZ = FALSE, EV = TRUE,
#' FA = FALSE, MASI = FALSE, MASV = TRUE, SIPC = TRUE,
#' ZA = FALSE)
ammistability <- function(model, n, alpha = 0.05,
ssi.method = c("farshadfar", "rao"), a = 1,
AMGE = TRUE, ASI = TRUE, ASV = TRUE, ASTAB = TRUE,
AVAMGE = TRUE, DA = TRUE, DZ = TRUE, EV = TRUE,
FA = TRUE, MASI = TRUE, MASV = TRUE, SIPC = TRUE,
ZA = TRUE, force.grouping = TRUE,
line.size = 1, line.alpha = 0.5, line.col = NULL,
point.size = 1, point.alpha = 0.5, point.col = NULL,
text.size = 2) {
# Check model class
if (!is(model, "AMMI")) {
stop('"model" is not of class "AMMI"')
}
# Check alpha value
if (!(0 < alpha && alpha < 1)) {
stop('"alpha" should be between 0 and 1 (0 < alpha < 1)')
}
# Find number of significant IPCs according to F test
if (missing(n) || is.null(n)) {
n <- sum(model$analysis$Pr.F <= alpha, na.rm = TRUE)
}
# Check for n
if (n %% 1 != 0 && length(n) != 1) {
stop('"n" is not an integer vector of unit length')
}
ssi.method <- match.arg(ssi.method)
sps <- data.frame(sp = c("AMGE", "ASI", "ASV", "ASTAB", "AVAMGE", "DA", "DZ",
"EV", "FA", "MASI", "MASV", "SIPC", "ZA"),
logical = c(AMGE, ASI, ASV, ASTAB, AVAMGE, DA, DZ,
EV, FA, MASI, MASV, SIPC, ZA))
if (nrow(sps[sps$logical, ]) < 2) {
stop('At least two stability parameters should be selected for computation')
}
means <- ASI.AMMI(model)
means$genotype <- rownames(means)
rownames(means) <- NULL
means <- means[, c("genotype", "means")]
SP <- means
SSI <- means
if (AMGE == TRUE) {
res <- AMGE.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "AMGE")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "AMGE"
rm(res)
}
if (ASI == TRUE) {
res <- MASI.AMMI(model, n = 2, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
colnames(res)[colnames(res) == "MASI"] <- "ASI"
SP <- merge(SP, res[, c("genotype", "ASI")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "ASI"
rm(res)
}
if (ASV == TRUE) {
res <- MASV.AMMI(model, n = 2, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
colnames(res)[colnames(res) == "MASV"] <- "ASV"
SP <- merge(SP, res[, c("genotype", "ASV")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "ASV"
rm(res)
}
if (ASTAB == TRUE) {
res <- ASTAB.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "ASTAB")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "ASTAB"
rm(res)
}
if (AVAMGE == TRUE) {
res <- AVAMGE.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "AVAMGE")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "AVAMGE"
rm(res)
}
if (DA == TRUE) {
res <- DA.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "DA")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "DA"
rm(res)
}
if (DZ == TRUE) {
res <- DZ.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "DZ")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "DZ"
rm(res)
}
if (EV == TRUE) {
res <- EV.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "EV")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "EV"
rm(res)
}
if (FA == TRUE) {
res <- FA.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "FA")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "FA"
rm(res)
}
if (MASI == TRUE) {
res <- MASI.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "MASI")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "MASI"
rm(res)
}
if (MASV == TRUE) {
res <- MASV.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "MASV")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "MASV"
rm(res)
}
if (SIPC == TRUE) {
res <- SIPC.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "SIPC")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "SIPC"
rm(res)
}
if (ZA == TRUE) {
res <- ZA.AMMI(model, n = n, alpha = alpha,
ssi.method = ssi.method, a = a)
res$genotype <- rownames(res)
SP <- merge(SP, res[, c("genotype", "Za")], by = "genotype")
SSI <- merge(SSI, res[, c("genotype", "SSI")], by = "genotype")
colnames(SSI)[colnames(SSI) == "SSI"] <- "Za"
rm(res)
}
# Correlations
#-----------------------------------------------------------------------------
# SP corr
del <- names(SP) %in% c("genotype", "means")
spcorr <- cor(SP[, !del])
sppmat <- ggcorrplot::cor_pmat(SP[, !del])
spcorrdf <- formatC(round(spcorr, 2), digits = 2, format = "f")
spcorrdf[] <- paste0(spcorrdf,
ifelse(sppmat < .01, "**", ifelse(sppmat < .05, "*", "")),
sep = "")
spcorrdf[upper.tri(spcorrdf)] <- NA
spcorrdf <- data.frame(spcorrdf, stringsAsFactors = FALSE)
spcorrg <- ggcorrplot(spcorr, hc.order = FALSE, type = "lower",
outline.color = "white", p.mat = sppmat,
ggtheme = theme_bw, show.diag = TRUE, lab = TRUE)
spcorrg <- spcorrg +
ggtitle("Correlation between different\nAMMI stability parameters")
spcorrg <- spcorrg +
theme(axis.text = element_text(colour = "black"))
# SSI corr
del <- names(SSI) %in% c("genotype", "means")
ssicorr <- cor(SSI[, !del])
ssipmat <- ggcorrplot::cor_pmat(SSI[, !del])
ssicorrdf <- formatC(round(ssicorr, 2), digits = 2, format = "f")
ssicorrdf[] <- paste0(ssicorrdf,
ifelse(ssipmat < .01, "**",
ifelse(ssipmat < .05, "*", "")),
sep = "")
ssicorrdf[upper.tri(ssicorrdf)] <- NA
ssicorrdf <- data.frame(ssicorrdf, stringsAsFactors = FALSE)
ssicorrg <- ggcorrplot(ssicorr, hc.order = FALSE, type = "lower",
outline.color = "white", p.mat = ssipmat,
ggtheme = theme_bw, show.diag = TRUE, lab = TRUE)
ssicorrg <- ssicorrg +
ggtitle("Correlation between simultaneous selection indices\nfrom different AMMI stability parameters")
ssicorrg <- ssicorrg +
theme(axis.text = element_text(colour = "black"))
# SP SSI corr
colnames(SSI) <- paste(colnames(SSI), "_SSI", sep = "")
colnames(SSI)[colnames(SSI) == "genotype_SSI"] <- "genotype"
colnames(SSI)[colnames(SSI) == "means_SSI"] <- "means"
spssi <- merge(SP, SSI[, names(SSI)[!(names(SSI) %in% "means")]],
by = "genotype")
del <- names(spssi) %in% c("genotype", "means")
spssicorr <- cor(spssi[, !del])
spssipmat <- ggcorrplot::cor_pmat(spssi[, !del])
spssicorrdf <- formatC(round(spssicorr, 2), digits = 2, format = "f")
spssicorrdf[] <- paste0(spssicorrdf,
ifelse(spssipmat < .01, "**",
ifelse(spssipmat < .05, "*", "")),
sep = "")
spssicorrdf[upper.tri(spssicorrdf)] <- NA
spssicorrdf <- data.frame(spssicorrdf, stringsAsFactors = FALSE)
spssicorrg <- ggcorrplot(spssicorr, hc.order = FALSE, type = "lower",
outline.color = "white", p.mat = spssipmat,
ggtheme = theme_bw, show.diag = TRUE, lab = TRUE,
title = )
spssicorrg <- spssicorrg +
ggtitle("Correlation between different AMMI stability parameters\nand corresponding simultaneous selection indices")
spssicorrg <- spssicorrg +
theme(axis.text = element_text(colour = "black"))
# Bump plots/ slopegraph
#-----------------------------------------------------------------------------
SPrank <- rankdf(df = SP, decreasing = "means",
increasing = names(SP)[!(names(SP) %in% c("genotype",
"means"))])
SPslopeg <- rankslopegraph(df = SPrank, names = "genotype",
force.grouping = TRUE, point.size = point.size,
point.alpha = point.alpha, line.size = line.size,
line.alpha = line.alpha, text.size = text.size,
legend.position = "none") +
ggtitle("Slopegraph of ranks of mean yields and AMMI stability parameters")
SPslopeg <- SPslopeg +
theme(axis.text = element_text(colour = "black"))
SSIrank <- rankdf(df = SSI, decreasing = "means",
increasing = names(SSI)[!(names(SSI) %in% c("genotype",
"means"))])
SSIslopeg <- rankslopegraph(df = SSIrank, names = "genotype",
force.grouping = TRUE, point.size = point.size,
point.alpha = point.alpha, line.size = line.size,
line.alpha = line.alpha, text.size = text.size,
legend.position = "none") +
ggtitle("Slopegraph of ranks of mean yields and simultaneous selction indices")
SSIslopeg <- SSIslopeg +
theme(axis.text = element_text(colour = "black"))
# heatmap
#-----------------------------------------------------------------------------
SPrank$genotype <- reorder(SPrank$genotype, SPrank$means)
SPrankmelt <- melt(SPrank[order(SPrank$means),], id.vars = "genotype")
SPheat <- ggplot(SPrankmelt, aes(x = variable, y = genotype)) +
geom_tile(aes(fill = value)) +
scale_x_discrete(expand = c(0, 0)) +
scale_y_discrete(expand = c(0, 0)) +
scale_fill_gradient(low = "black", high = "green") +
ylab("Genotype") + labs(fill = "Rank") +
theme_bw() +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
axis.text = element_text(colour = "black"),
axis.title.x = element_blank(),
legend.position = "bottom")
SSIrank$genotype <- reorder(SSIrank$genotype, SSIrank$means)
SSIrankmelt <- melt(SSIrank[order(SSIrank$means),], id.vars = "genotype")
SSIheat <- ggplot(SSIrankmelt, aes(x = variable, y = genotype)) +
geom_tile(aes(fill = value)) +
scale_x_discrete(expand = c(0, 0)) +
scale_y_discrete(expand = c(0, 0)) +
scale_fill_gradient(low = "black", high = "green") +
ylab("Genotype") + labs(fill = "Rank") +
theme_bw() +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
axis.text = element_text(colour = "black"),
axis.title.x = element_blank(),
legend.position = "bottom")
details <- list(`Stability parameters estimated` = as.character(sps[sps$logical, ]$sp),
`SSI method` = ifelse(ssi.method == "rao",
paste("Rao and Prabhakaran (2005)",
"with a =", a),
"Farshadfar (2008)"))
out <- list(Details = details,
`Stability Parameters` = SP,
`Simultaneous Selection Indices` = SSI,
`SP Correlation` = spcorrdf,
`SSI Correlation` = ssicorrdf,
`SP and SSI Correlation` = spssicorrdf,
`SP Correlogram` = spcorrg,
`SSI Correlogram` = ssicorrg,
`SP and SSI Correlogram` = spssicorrg,
`SP Slopegraph` = SPslopeg,
`SSI Slopegraph` = SSIslopeg,
`SP Heatmap` = SPheat,
`SSI Heatmap` = SSIheat)
return(out)
}
Try the ammistability package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.