jranke / mkin
1
#' Plot model fits (observed and fitted) and the residuals for a row or column
2
#' of an mmkin object
3
#'
4
#' When x is a row selected from an mmkin object (\code{\link{[.mmkin}}), the
5
#' same model fitted for at least one dataset is shown. When it is a column,
6
#' the fit of at least one model to the same dataset is shown.
7
#'
8
#' If the current plot device is a \code{\link[tikzDevice]{tikz}} device, then
9
#' latex is being used for the formatting of the chi2 error level.
10
#'
11
#' @param x An object of class \code{\link{mmkin}}, with either one row or one
12
#'   column.
13
#' @param main The main title placed on the outer margin of the plot.
14
#' @param legends An index for the fits for which legends should be shown.
15
#' @param resplot Should the residuals plotted against time, using
16
#'   \code{\link{mkinresplot}}, or as squared residuals against predicted
17
#'   values, with the error model, using \code{\link{mkinerrplot}}.
18
#' @param ylab Label for the y axis.
19
#' @param standardized Should the residuals be standardized? This option
20
#'   is passed to \code{\link{mkinresplot}}, it only takes effect if
21
#'   `resplot = "time"`.
22
#' @param show_errmin Should the chi2 error level be shown on top of the plots
23
#'   to the left?
24
#' @param errmin_var The variable for which the FOCUS chi2 error value should
25
#'   be shown.
26
#' @param errmin_digits The number of significant digits for rounding the FOCUS
27
#'   chi2 error percentage.
28
#' @param cex Passed to the plot functions and \code{\link{mtext}}.
29
#' @param rel.height.middle The relative height of the middle plot, if more
30
#'   than two rows of plots are shown.
31
#' @param ymax Maximum y axis value for \code{\link{plot.mkinfit}}.
32
#' @param \dots Further arguments passed to \code{\link{plot.mkinfit}} and
33
#'   \code{\link{mkinresplot}}.
34
#' @return The function is called for its side effect.
35
#' @author Johannes Ranke
36
#' @examples
37
#'
38
#'   \dontrun{
39
#'   # Only use one core not to offend CRAN checks
40
#'   fits <- mmkin(c("FOMC", "HS"),
41
#'                 list("FOCUS B" = FOCUS_2006_B, "FOCUS C" = FOCUS_2006_C), # named list for titles
42
#'                 cores = 1, quiet = TRUE, error_model = "tc")
43
#'   plot(fits[, "FOCUS C"])
44
#'   plot(fits["FOMC", ])
45
#'   plot(fits["FOMC", ], show_errmin = FALSE)
46
#'
47
#'   # We can also plot a single fit, if we like the way plot.mmkin works, but then the plot
48
#'   # height should be smaller than the plot width (this is not possible for the html pages
49
#'   # generated by pkgdown, as far as I know).
50
#'   plot(fits["FOMC", "FOCUS C"]) # same as plot(fits[1, 2])
51
#'
52
#'   # Show the error models
53
#'   plot(fits["FOMC", ], resplot = "errmod")
54
#'   }
55
#'
56
#' @export
57
plot.mmkin <- function(x, main = "auto", legends = 1,
58
  resplot = c("time", "errmod"),
59
  ylab = "Residue",
60
  standardized = FALSE,
61
  show_errmin = TRUE,
62
  errmin_var = "All data", errmin_digits = 3,
63
  cex = 0.7, rel.height.middle = 0.9,
64
  ymax = "auto", ...)
65
{
66

67 1
  oldpar <- par(no.readonly = TRUE)
68 1
  on.exit(par(oldpar, no.readonly = TRUE))
69

70 1
  n.m <- nrow(x)
71 1
  n.d <- ncol(x)
72

73 1
  resplot <- match.arg(resplot)
74

75
  # We can handle either a row (different models, same dataset)
76
  # or a column (same model, different datasets)
77 0
  if (n.m > 1 & n.d > 1) stop("Please select fits either for one model or for one dataset")
78 0
  if (n.m == 1 & n.d == 1) loop_over = "none"
79 1
  if (n.m > 1) loop_over <- "models"
80 1
  if (n.d > 1) loop_over <- "datasets"
81 1
  n.fits <- length(x)
82

83
  # Set the main plot titles from the names of the models or the datasets
84
  # Will be integer indexes if no other names are present in the mmkin object
85 1
  if (main == "auto") {
86 1
    main = switch(loop_over,
87 1
                  none = paste(rownames(x), colnames(x)),
88 1
                  models = colnames(x),
89 1
                  datasets = rownames(x))
90
  }
91

92
  # Set relative plot heights, so the first and the last plot are the norm
93
  # and the middle plots (if n.fits >2) are smaller by rel.height.middle
94 1
  rel.heights <- if (n.fits > 2) c(1, rep(rel.height.middle, n.fits - 2), 1)
95 1
                 else rep(1, n.fits)
96 1
  layout(matrix(1:(2 * n.fits), n.fits, 2, byrow = TRUE), heights = rel.heights)
97

98 1
  par(cex = cex)
99

100 1
  for (i.fit in 1:n.fits) {
101

102
    # Margins for top row of plots when we have more than one row
103
    # Reduce bottom margin by 2.1 - hides x axis legend
104 1
    if (i.fit == 1 & n.fits > 1) {
105 1
      par(mar = c(3.0, 4.1, 4.1, 2.1))
106
    }
107

108
    # Margins for middle rows of plots, if any
109 1
    if (i.fit > 1 & i.fit < n.fits) {
110
      # Reduce top margin by 2 after the first plot as we have no main title,
111
      # reduced plot height, therefore we need rel.height.middle in the layout
112 1
      par(mar = c(3.0, 4.1, 2.1, 2.1))
113
    }
114

115
    # Margins for bottom row of plots when we have more than one row
116 1
    if (i.fit == n.fits & n.fits > 1) {
117
      # Restore bottom margin for last plot to show x axis legend
118 1
      par(mar = c(5.1, 4.1, 2.1, 2.1))
119
    }
120

121 1
    fit <- x[[i.fit]]
122 1
    if (ymax == "auto") {
123 1
      plot(fit, legend = legends == i.fit, ylab = ylab, ...)
124
    } else {
125 0
      plot(fit, legend = legends == i.fit, ylim = c(0, ymax), ylab = ylab, ...)
126
    }
127

128 1
    title(main, outer = TRUE, line = -2)
129

130 1
    fit_name <- switch(loop_over,
131 1
                       models = rownames(x)[i.fit],
132 1
                       datasets = colnames(x)[i.fit],
133 1
                       none = "")
134

135 1
    if (show_errmin) {
136 1
      chi2 <- signif(100 * mkinerrmin(fit)[errmin_var, "err.min"], errmin_digits)
137

138
      # Use LateX if the current plotting device is tikz
139 1
      if (names(dev.cur()) == "tikz output") {
140 0
        chi2_text <- paste0(fit_name, " $\\chi^2$ error level = ", chi2, "\\%")
141
      } else {
142 1
        chi2_perc <- paste0(chi2, "%")
143 1
        chi2_text <- bquote(.(fit_name) ~ chi^2 ~ "error level" == .(chi2_perc))
144
      }
145 1
      mtext(chi2_text, cex = cex, line = 0.4)
146
    } else {
147 0
      mtext(fit_name, cex = cex, line = 0.4)
148
    }
149

150 1
    if (resplot == "time") {
151 1
      mkinresplot(fit, legend = FALSE, standardized = standardized, ...)
152
    } else {
153 0
      mkinerrplot(fit, legend = FALSE, ...)
154
    }
155 1
    mtext(paste(fit_name, "residuals"), cex = cex, line = 0.4)
156
  }
157
}

Read our documentation on viewing source code .

Loading