#' Plot model fits (observed and fitted) and the residuals for a row or column
#' of an mmkin object
#'
#' When x is a row selected from an mmkin object (\code{\link{[.mmkin}}), the
#' same model fitted for at least one dataset is shown. When it is a column,
#' the fit of at least one model to the same dataset is shown.
#'
#' If the current plot device is a \code{\link[tikzDevice]{tikz}} device, then
#' latex is being used for the formatting of the chi2 error level.
#'
#' @param x An object of class \code{\link{mmkin}}, with either one row or one
#' column.
#' @param main The main title placed on the outer margin of the plot.
#' @param legends An index for the fits for which legends should be shown.
#' @param resplot Should the residuals plotted against time, using
#' \code{\link{mkinresplot}}, or as squared residuals against predicted
#' values, with the error model, using \code{\link{mkinerrplot}}.
#' @param show_errmin Should the chi2 error level be shown on top of the plots
#' to the left?
#' @param errmin_var The variable for which the FOCUS chi2 error value should
#' be shown.
#' @param errmin_digits The number of significant digits for rounding the FOCUS
#' chi2 error percentage.
#' @param cex Passed to the plot functions and \code{\link{mtext}}.
#' @param rel.height.middle The relative height of the middle plot, if more
#' than two rows of plots are shown.
#' @param \dots Further arguments passed to \code{\link{plot.mkinfit}} and
#' \code{\link{mkinresplot}}.
#' @return The function is called for its side effect.
#' @author Johannes Ranke
#' @examples
#'
#' \dontrun{
#' # Only use one core not to offend CRAN checks
#' fits <- mmkin(c("FOMC", "HS"),
#' list("FOCUS B" = FOCUS_2006_B, "FOCUS C" = FOCUS_2006_C), # named list for titles
#' cores = 1, quiet = TRUE, error_model = "tc")
#' plot(fits[, "FOCUS C"])
#' plot(fits["FOMC", ])
#'
#' # We can also plot a single fit, if we like the way plot.mmkin works, but then the plot
#' # height should be smaller than the plot width (this is not possible for the html pages
#' # generated by pkgdown, as far as I know).
#' plot(fits["FOMC", "FOCUS C"]) # same as plot(fits[1, 2])
#'
#' # Show the error models
#' plot(fits["FOMC", ], resplot = "errmod")
#' }
#'
#' @export
plot.mmkin <- function(x, main = "auto", legends = 1,
resplot = c("time", "errmod"),
show_errmin = TRUE,
errmin_var = "All data", errmin_digits = 3,
cex = 0.7, rel.height.middle = 0.9, ...) {
n.m <- nrow(x)
n.d <- ncol(x)
resplot <- match.arg(resplot)
# We can handle either a row (different models, same dataset)
# or a column (same model, different datasets)
if (n.m > 1 & n.d > 1) stop("Please select fits either for one model or for one dataset")
if (n.m == 1 & n.d == 1) loop_over = "none"
if (n.m > 1) loop_over <- "models"
if (n.d > 1) loop_over <- "datasets"
n.fits <- length(x)
# Set the main plot titles from the names of the models or the datasets
# Will be integer indexes if no other names are present in the mmkin object
if (main == "auto") {
main = switch(loop_over,
none = paste(rownames(x), colnames(x)),
models = colnames(x),
datasets = rownames(x))
}
oldpar <- par(no.readonly = TRUE)
# Set relative plot heights, so the first and the last plot are the norm
# and the middle plots (if n.fits >2) are smaller by rel.height.middle
rel.heights <- if (n.fits > 2) c(1, rep(rel.height.middle, n.fits - 2), 1)
else rep(1, n.fits)
layout(matrix(1:(2 * n.fits), n.fits, 2, byrow = TRUE), heights = rel.heights)
par(cex = cex)
for (i.fit in 1:n.fits) {
# Margins for top row of plots when we have more than one row
# Reduce bottom margin by 2.1 - hides x axis legend
if (i.fit == 1 & n.fits > 1) {
par(mar = c(3.0, 4.1, 4.1, 2.1))
}
# Margins for middle rows of plots, if any
if (i.fit > 1 & i.fit < n.fits) {
# Reduce top margin by 2 after the first plot as we have no main title,
# reduced plot height, therefore we need rel.height.middle in the layout
par(mar = c(3.0, 4.1, 2.1, 2.1))
}
# Margins for bottom row of plots when we have more than one row
if (i.fit == n.fits & n.fits > 1) {
# Restore bottom margin for last plot to show x axis legend
par(mar = c(5.1, 4.1, 2.1, 2.1))
}
fit <- x[[i.fit]]
plot(fit, legend = legends == i.fit, ...)
title(main, outer = TRUE, line = -2)
fit_name <- switch(loop_over,
models = rownames(x)[i.fit],
datasets = colnames(x)[i.fit],
none = "")
if (show_errmin) {
chi2 <- signif(100 * mkinerrmin(fit)[errmin_var, "err.min"], errmin_digits)
# Use LateX if the current plotting device is tikz
if (names(dev.cur()) == "tikz output") {
chi2_text <- paste0(fit_name, " $\\chi^2$ error level = ", chi2, "\\%")
} else {
chi2_perc <- paste0(chi2, "%")
chi2_text <- bquote(.(fit_name) ~ chi^2 ~ "error level" == .(chi2_perc))
}
mtext(chi2_text, cex = cex, line = 0.4)
}
if (resplot == "time") {
mkinresplot(fit, legend = FALSE, ...)
} else {
mkinerrplot(fit, legend = FALSE, ...)
}
if (show_errmin) {
mtext(paste(fit_name, "residuals"), cex = cex, line = 0.4)
}
}
par(oldpar, no.readonly = TRUE)
}
