/
plot.mmkin.R
157 lines (140 loc) · 5.92 KB
/
plot.mmkin.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
#' 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 ylab Label for the y axis.
#' @param standardized Should the residuals be standardized? This option
#' is passed to \code{\link{mkinresplot}}, it only takes effect if
#' `resplot = "time"`.
#' @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 ymax Maximum y axis value for \code{\link{plot.mkinfit}}.
#' @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", ])
#' plot(fits["FOMC", ], show_errmin = FALSE)
#'
#' # 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"),
ylab = "Residue",
standardized = FALSE,
show_errmin = TRUE,
errmin_var = "All data", errmin_digits = 3,
cex = 0.7, rel.height.middle = 0.9,
ymax = "auto", ...)
{
oldpar <- par(no.readonly = TRUE)
on.exit(par(oldpar, no.readonly = TRUE))
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))
}
# 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]]
if (ymax == "auto") {
plot(fit, legend = legends == i.fit, ylab = ylab, ...)
} else {
plot(fit, legend = legends == i.fit, ylim = c(0, ymax), ylab = ylab, ...)
}
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)
} else {
mtext(fit_name, cex = cex, line = 0.4)
}
if (resplot == "time") {
mkinresplot(fit, legend = FALSE, standardized = standardized, ...)
} else {
mkinerrplot(fit, legend = FALSE, ...)
}
mtext(paste(fit_name, "residuals"), cex = cex, line = 0.4)
}
}