Nothing
#' @title
#' Plotting for microsynth objects.
#'
#' @description
#' Using a \code{microsynth} object as an input, this function gives time
#' series plots of selected outcomes.
#'
#' @details
#' Plots are given over both pre- and
#' intervention time periods and shown in terms of raw outcome values or
#' treatment/control differences. Time series of permutation groups may be
#' overlaid to help illustrate statistical uncertainty.
#'
#' Only required input is a parameter \code{ms} which is a microsynth object.
#'
#' @param ms A microsynth object
#'
#' @param plot.var A vector of variable names giving the outcome variables that
#' are shown in plots. If \code{plot.var = NULL}, all outcome variables that
#' are included in \code{ms} are plotted. Only variables contained in
#' the input \code{result.var} as used in the creation of \code{ms} can
#' be plotted using \code{plot()}.
#'
#' @param file A character string giving the name of file that will be
#' created in the home directory containing plots.
#' The name should have a \code{.pdf} extension.
#'
#' @param start.pre An integer indicating the time point that corresponds to the
#' earliest time period that will be plotted.
#' When \code{start.pre = NULL}, it is reset to the
#' minimum time appearing in \code{ms}.
#'
#' @param end.pre An integer that gives the final time point of the
#' pre-intervention period. That is, \code{end.pre} is the last time at
#' which treatment and synthetic control will were matched to one another.
#' All time points
#' following \code{end.pre} are considered to be post-intervention and the
#' behavior of outcomes will be compared between the treatment and synthetic
#' control groups across those time periods.
#' If \code{end.pre = NULL} the end of the pre-intervention period will be
#' determined from the object \code{ms}.
#'
#' @param end.post An integer that gives final time point that will be plotted.
#' When \code{end.post = NULL} (the default), it is reset
#' to the maximum time that appears in \code{ms}.
#'
#' @param sep If \code{sep = TRUE}, separate plots will be generated for each
#' outcome. Applicable only if plots are saved to file (
#' \code{plot.file} is \code{non-NULL}). To change display of plots produced
#' as output, use \code{\link[graphics]{par}}.
#'
#' @param legend.spot The location of the legend in the plots.
#'
#' @param plot.first The number of permutation groups to plot.
#'
#' @param height The height of the graphics region (in inches)
#' when a pdf is created.
#'
#' @param width The width of the graphics region (in inches)
#' when a pdf is created.
#'
#' @param at A vector that gives the location of user-specified x-axis labels.
#' \code{at} should be a (numeric) subset of the named time points contained
#' in \code{ms} (e.g., \code{colnames(ms$Plot.Stats$Treatment)}).
#'
#' @param labels A vector of the same length as \code{at} that gives the names
#' of the labels that will be marked at the times indicated by \code{at} in
#' the plots.
#'
#' @param all A scalar character string giving the unit name for cases.
#' If \code{NULL}, a third curve showing the overall outcome levels is
#' not plotted.
#'
#' @param main.tc A scalar (or a vector of the same length as \code{plot.var})
#' character string giving the title to be used for the first plots
#' (that show treatment and control). Defaults to \code{plot.var}.
#'
#' @param main.diff A scalar (or a vector of the same length as \code{plot.var})
#' character string giving the title to be used for the second plots
#' (that show differences between treatment and control).
#' Defaults to \code{plot.var}.
#'
#' @param xlab.tc A scalar (or a vector of the same length as \code{plot.var})
#' character string giving the x-axis labels to be used for the first plots
#' (that show treatment and control). Defaults to \code{''}.
#'
#' @param xlab.diff A scalar (or a vector of the same length as \code{plot.var})
#' character string giving the x-axis labels to be used for the second plots
#' (that show differences between treatment and control).
#' Defaults to \code{''}.
#'
#' @param ylab.tc A scalar (or a vector of the same length as \code{plot.var})
#' character string giving the y-axis labels to be used for the first plots
#' (that show treatment and control). Defaults to \code{plot.var}.
#'
#' @param ylab.diff A scalar (or a vector of the same length as \code{plot.var})
#' character string giving the y-axis labels to be used for the second plots
#' (that show differences between treatment and control).
#' Defaults to \code{'Treatment - Control'}.
#'
#' @return No return value, called for side effects (i.e., to produce plots
#' of outcome values and treatment/control differences, with the option to
#' write to file).
#'
#' @examples
#'
#' # Declare time-variant (outcome) and time-invariant variables for matching
#' cov.var <- c('TotalPop', 'BLACK', 'HISPANIC', 'Males_1521',
#' 'HOUSEHOLDS', 'FAMILYHOUS', 'FEMALE_HOU', 'RENTER_HOU', 'VACANT_HOU')
#'
#' match.out <- c('i_felony', 'i_misdemea', 'i_drugs', 'any_crime')
#'
#' set.seed(99199) # for reproducibility
#'
#' \donttest{
#'
#' # Perform matching and estimation, without permutations or jackknife
#' # runtime: <1 min
#' sea1 <- microsynth(seattledmi,
#' idvar='ID', timevar='time', intvar='Intervention',
#' start.pre=1, end.pre=12, end.post=16,
#' match.out=match.out, match.covar=cov.var,
#' result.var=match.out, omnibus.var=match.out,
#' test='lower',
#' n.cores = min(parallel::detectCores(), 2))
#'
#' # Plot with default settings in the GUI.
#' plot_microsynth(sea1)
#'
#' # Make plots, display, and save to a single file (plots.pdf).
#' plot_microsynth(sea1, file = file.path(tempdir(), 'plots.pdf'), sep = FALSE)
#'
#' # Make plots for only one outcome, display, and save to a single file.
#' plot_microsynth(sea1, plot.var = 'any_crime',
#' file = file.path(tempdir(), 'plots.pdf'), sep = FALSE)
#' }
#'
#' @export
plot_microsynth <- function(ms, plot.var = NULL, start.pre = NULL, end.pre = NULL, end.post = NULL, file = NULL, sep = TRUE, plot.first = NULL, legend.spot = "bottomleft", height = NULL, width = NULL, at = NULL, labels = NULL,
all = "cases", main.tc = NULL, main.diff = NULL, xlab.tc = NULL, xlab.diff = NULL, ylab.tc = NULL, ylab.diff = NULL) {
if (!is.element("Plot.Stats", names(ms))) {
stop("object ms does not contain output regarding results (e.g., only weights were generated).")
}
plotdat.t <- ms$Plot.Stats$Treatment
plotdat.c <- ms$Plot.Stats$Control
if (length(all) > 0) {
plotdat.a <- ms$Plot.Stats$All
} else {
plotdat.a <- ms$Plot.Stats$Treatment
}
plotdat.d <- ms$Plot.Stats$Difference
scale.by <- ms$Plot.Stats$scale.by
if (length(height) == 0) {
if (!sep) {
height <- 11
} else {
height <- 5
}
}
if (length(width) == 0) {
if (!sep) {
width <- 8.5
} else {
width <- 5
}
}
if (length(plot.first) == 0) {
plot.first <- dim(plotdat.d)[2] - 1
}
time.names <- colnames(plotdat.t)
if (length(start.pre) > 0) {
if (is.element(as.character(start.pre), time.names)) {
start.pre <- match(as.character(start.pre), time.names)
} else {
start.pre <- which(as.numeric(time.names) >= as.numeric(start.pre))[1]
}
} else {
start.pre <- 1
}
if (length(end.pre) > 0) {
if (is.element(as.character(end.pre), time.names)) {
end.pre <- match(as.character(end.pre), time.names)
} else {
end.pre <- which(as.numeric(time.names) <= as.numeric(end.pre))
end.pre <- end.pre[length(end.pre)]
}
} else {
end.pre <- ms$Plot.Stats$end.pre
end.pre <- match(as.character(end.pre), time.names)
}
if (length(end.post) > 0) {
if (is.element(as.character(end.post), time.names)) {
end.post <- match(as.character(end.post), time.names)
} else {
end.post <- which(as.numeric(time.names) <= as.numeric(end.post))
end.post <- end.post[length(end.post)]
}
} else {
end.post <- length(time.names)
}
if (start.pre == end.post) {
stop("Do not plot: start.pre == end.post")
}
all.vars <- rownames(plotdat.t)
if (length(plot.var) == 0) {
plot.var <- all.vars
}
plot.var <- remove.vars(plot.var, all.vars, "plot.var")
if (length(xlab.tc) == 0) {
xlab.tc <- rep("", length(plot.var))
} else if (length(xlab.tc) == 1) {
xlab.tc <- rep(xlab.tc, length(plot.var))
}
if (length(xlab.diff) == 0) {
xlab.diff <- rep("", length(plot.var))
} else if (length(xlab.diff) == 1) {
xlab.diff <- rep(xlab.diff, length(plot.var))
}
if (length(ylab.tc) == 0) {
ylab.tc <- plot.var
} else if (length(ylab.tc) == 1) {
ylab.tc <- rep(ylab.tc, length(plot.var))
}
if (length(ylab.diff) == 0) {
ylab.diff <- rep("Treatment - Control", length(plot.var))
} else if (length(ylab.diff) == 1) {
ylab.diff <- rep(ylab.diff, length(plot.var))
}
if (length(main.tc) == 0) {
main.tc <- plot.var
} else if (length(main.tc) == 1) {
main.tc <- rep(main.tc, length(plot.var))
}
if (length(main.diff) == 0) {
main.diff <- plot.var
} else if (length(main.diff) == 1) {
main.diff <- rep(main.diff, length(plot.var))
}
# xnams <- as.numeric(colnames(test1))
xnams <- 1:length(time.names)
tuse <- xnams <= end.post & xnams >= start.pre
use <- xnams <= end.pre
nuse <- xnams >= end.pre
if (end.post > end.pre) {
fuse <- !use & xnams <= end.post
} else if (end.post == end.pre) {
fuse <- xnams >= end.pre & xnams <= end.post
} else {
stop("end.post is less than end.pre")
}
oldpar <- graphics::par(no.readonly = TRUE)
on.exit(graphics::par(oldpar))
if (length(file) == 0) {
graphics::par(mfrow = c(1, 2), ask = FALSE)
} else {
file.type <- "pdf"
if (substr(file, nchar(file) - 3, nchar(file)) == ".pdf") {
file <- substr(file, 1, nchar(file) - 4)
} else if (substr(file, nchar(file) - 3, nchar(file)) == ".png") {
file <- substr(file, 1, nchar(file) - 4)
file.type <- "png"
}
if (!sep) {
if (file.type == "pdf") {
file <- paste(file, ".pdf", sep = "")
grDevices::pdf(file = file, width = width, height = height)
} else if (file.type == "png") {
file <- paste(file, ".png", sep = "")
grDevices::png(file = file, width = width, height = height, units = "in", res = 500)
}
graphics::par(mfrow = c(3, 2), ask = FALSE)
}
}
for (j in 1:length(plot.var)) {
# use.mu <- which(mu[, plot.var[j]] > 0)
use.mu <- 1:dim(plotdat.d)[2]
for (i in 1:dim(plotdat.d)[2]) {
tmp <- plotdat.d[plot.var[j], i, ]
if (i == 1) {
if (sep & length(file) > 0) {
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, "_", plot.var[j], "_TC.pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, "_", plot.var[j], "_TC.png", sep = ""), width = width, height = height, units = "in", res = 500)
}
}
tmp1 <- plotdat.t[plot.var[j], ]
tmp2 <- plotdat.c[plot.var[j], ]
tmp3 <- scale.by * plotdat.a[plot.var[j], ]
lty <- c(1, 2, 4)
col <- c(2, 1, 3)
lwd <- c(2, 2, 2)
ylim <- c(min(tmp1[tuse], tmp2[tuse]), max(tmp1[tuse], tmp2[tuse]))
ylim <- c(min(ylim, tmp3[tuse]), max(ylim, tmp3[tuse]))
ylim[2] <- 1.2 * ylim[2]
xxnams1 <- as.numeric(as.character(time.names[xnams]))
iend.pre <- as.numeric(as.character(time.names[end.pre]))
if (sum(is.na(xxnams1)) > 0) {
xxnams1 <- xnams
iend.pre <- end.pre
}
if (length(at) == 0) {
# at <- xxnams1
xaxt <- "s"
} else {
xaxt <- "n"
}
# if(length(labels) == 0) { labels <- at }
xxnams2 <- xxnams1[tuse & use]
xxnams3 <- xxnams1[tuse & nuse]
xxnams1 <- xxnams1[tuse]
xlim <- c(min(xxnams1), max(xxnams1))
graphics::plot(xxnams1, tmp1[tuse], type = "l", lty = lty[1], col = col[1], lwd = lwd[1], xlim = xlim, xlab = xlab.tc[j], ylab = ylab.tc[j], main = main.tc[j], ylim = ylim, xaxt = xaxt)
if (length(at) > 0) {
at1 <- intersect(at, xxnams1)
labels1 <- labels[is.element(at, xxnams1)]
graphics::axis(1, at = at1, labels = labels1)
}
graphics::abline(v = iend.pre, lty = 2, col = 2)
graphics::lines(xxnams1, tmp2[tuse], type = "l", lty = lty[2], col = col[2], lwd = lwd[2])
if (length(all) > 0) {
graphics::lines(xxnams1, tmp3[tuse], type = "l", lty = lty[3], col = col[3], lwd = lwd[3])
leg3 <- paste("All ", all, " (scaled)", sep = "")
leg <- c("Treatment", "Synthetic Control", leg3)
} else {
leg <- c("Treatment", "Synthetic Control")
col <- col[1:2]
lty <- lty[1:2]
lwd <- lwd[1:2]
}
graphics::legend(legend.spot, legend = leg, col = col, lty = lty, lwd = lwd, cex = 0.8, bty = "n")
if (sep & length(file) > 0) {
grDevices::dev.off()
}
if (is.element(i, use.mu)) {
if (sep & length(file) > 0) {
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, "_", plot.var[j], "_Diff.pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, "_", plot.var[j], "_Diff.png", sep = ""), width = width, height = height, units = "in", res = 500)
}
}
ylim1 <- c(min(tmp[tuse]), max(tmp[tuse]))
bigtmp <- plotdat.d[plot.var[j], use.mu, ]
bigtmp <- matrix(bigtmp, length(use.mu), dim(plotdat.d)[3])
ylim2a <- 2 * min(apply(bigtmp[, tuse, drop = FALSE], 2, stats::quantile, probs = 0.05, na.rm = TRUE))
ylim2b <- 2 * max(apply(bigtmp[, tuse, drop = FALSE], 2, stats::quantile, probs = 0.95, na.rm = TRUE))
ylim <- c(min(ylim1[1], ylim2a, na.rm = TRUE), max(ylim1[2], ylim2b, na.rm = TRUE))
xlim <- c(min(xxnams1), max(xxnams1))
graphics::plot(xxnams2, tmp[tuse & use], type = "l", ylim = ylim, xlim = xlim, col = 2, lty = 2, main = main.diff[j], xlab = xlab.diff[j], ylab = ylab.diff[j], xaxt = xaxt)
if (length(at) > 0) {
at1 <- intersect(at, xxnams1)
labels1 <- labels[is.element(at, xxnams1)]
graphics::axis(1, at = at1, labels = labels1)
}
} else {
if (!sep & length(file) > 0) {
graphics::plot(1, 1, main = main.diff[j])
}
}
} else {
if (i <= min((dim(plotdat.d)[2]), plot.first)) {
if (is.element(i, use.mu) & is.element(1, use.mu)) {
graphics::lines(xxnams1, tmp[tuse], col = "azure3")
}
}
}
if (i == dim(plotdat.d)[2] & is.element(1, use.mu)) {
tmp <- plotdat.d[plot.var[j], 1, ]
graphics::lines(xxnams2, tmp[tuse & use], lty = 1, col = 1, lwd = 2)
graphics::lines(xxnams3, tmp[tuse & nuse], col = 2, lwd = 2)
graphics::abline(v = iend.pre, col = 2, lwd = 1, lty = 2)
graphics::abline(h = 0, lwd = 1, lty = 2)
if (sep & length(file) > 0) {
grDevices::dev.off()
}
}
}
}
if (!sep & length(file) > 0) {
grDevices::dev.off()
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.