Nothing
plot_gerbil <- function(gerb, type = "Univariate", vars = NULL, file = NULL, sep = FALSE, imp = NULL,
obs.col = "blue4", imp.col = "brown2", obs.lty = 1, imp.lty = 2, obs.lwd = NULL,
imp.lwd = NULL, height = NULL, width = NULL, legend.spot = "topright", legend.text = legend.text,
mfrow = c(3, 2), trace.type = "Mean", log = NULL, pch = c(1, 3), partial = "imputed", ...)
{
if (tolower(type) == "univariate" | type == 1) {
plot_gerbil.uni(gerb = gerb, file = file, vars = vars, imp = imp, obs.col = obs.col, imp.col = imp.col, sep = sep,
obs.lty = obs.lty, imp.lty = imp.lty, use.log = log, legend.text = legend.text, legend.spot = legend.spot,
height = height, width = width, mfrow = mfrow, obs.lwd = obs.lwd, imp.lwd = imp.lwd, ...)
} else if (tolower(type) == "bivariate" | tolower(type) == "multivariate" | type == 2) {
#if (length(imp) > 1) {
# warning("imp can not have length greater than 1 for bivariate plotting. Will use only the first element.")
#}
plot_gerbil.mult(gerb = gerb, file = file, vars = vars, imp = imp, obs.col = obs.col, imp.col = imp.col,
lty = c(obs.lty, imp.lty), lwd = c(obs.lwd, imp.lwd), sep = sep, legend.text = legend.text,
legend.spot = legend.spot, height = height, width = width, mfrow = mfrow, pch = pch, use.log = log, partial = partial, ...)
} else if (tolower(type) == "ts" | tolower(type) == "time.series" | type == 3) {
plot_gerbil.ts(gerb = gerb, file = file, vars = vars, imps = imp, obs.col = obs.col, imp.col = imp.col, sep = sep, legend.text = legend.text,
legend.spot = legend.spot, height = height, width = width, mfrow = mfrow, type = trace.type,
obs.lty = obs.lty, imp.lty = imp.lty, obs.lwd = obs.lwd, imp.lwd = imp.lwd, ...)
} else {
stop("Input 'type' not recognized.")
}
}
plot_gerbil.ts <- function (gerb, file = NULL, vars = NULL, obs.col = "blue4", imp.col = "brown2", obs.lty = 2, imp.lty = 1,
imps = NULL, sep = FALSE, height = NULL, width = NULL,
type = "Mean", plot.obs = TRUE, obs.lwd = NULL, imp.lwd = NULL, mfrow = c(3, 2), legend.text = NULL,
legend.spot = "topright", main = NULL, xlab = NULL, ylab = NULL, bty = NULL, ...) {
oldpar <- par(no.readonly = TRUE)
on.exit(par(oldpar))
plot.legend <- TRUE
if (is.logical(legend.text)) {
if(!legend) {
plot.legend <- FALSE
} else {
legend.text <- NULL
}
} else if (length(legend.text) == 1) {
if (legend.text == "n") {
plot.legend <- FALSE
} else if (legend.text == "y") {
legend.text <- NULL
}
}
if (length(gerb$chainSeq) == 0) {
stop("The gerbil object was created with trace = FALSE. Cannot create trace plots.")
}
if (type == 1 | tolower(type) == "mean") {
gerb.tmp <- gerb$chainSeq$means.mis
gerb.tmp.obs <- gerb$chainSeq$means.obs
} else if (type == 2 | tolower(type) == "variance") {
if (is.element("vars.mis", names(gerb$chainSeq))) {
gerb.tmp <- gerb$chainSeq$vars.mis
gerb.tmp.obs <- gerb$chainSeq$vars.obs
} else {
stop("Input 'gerb' was produced with 'calc.var = FALSE'.")
}
} else if (type == 3 | tolower(type) == "latent.mean") {
gerb.tmp <- gerb$chainSeq$means.mis.latent
gerb.tmp.obs <- NULL
if (plot.obs) {
warning("Will not plot means of the latent process for observed values.")
plot.obs <- FALSE
}
} else if (type == 4 | tolower(type) == "latent.variance") {
if (is.element("vars.mis.latent", names(gerb$chainSeq))) {
gerb.tmp <- gerb$chainSeq$vars.mis.latent
gerb.tmp.obs <- NULL
if (plot.obs) {
warning("Will not plot variances of the latent process for observed values.")
plot.obs <- FALSE
}
} else {
stop("Input 'gerb' was produced with 'calc.var = FALSE'.")
}
} else {
stop ("Invalid value of input 'type'.")
}
if (length(vars) == 0) {
#vars <- dimnames(gerb.tmp)[[1]]
vars <- colnames(gerb$missing.latent)[colSums(gerb$missing.latent == 1 | gerb$missing.latent == 2) > 0]
}
if (length(file) == 0 | length(vars) == 1) {
sep <- TRUE
supp.suff <- TRUE
} else {
supp.suff <- FALSE
}
if (sep) {
mfrow <- c(1, 1)
}
if(length(height) == 0) {
if (!sep) {
height <- 11
} else {
height <- 6.5
}
}
if(length(width) == 0) {
if (!sep) {
width <- 8.5
} else {
width <- 6.5
}
}
if (length(file) == 0) {
graphics::par(mfrow = c(1, 1), 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 = mfrow, ask = FALSE)
}
}
latent.vars <- dimnames(gerb.tmp)[[1]]
orig.vars <- dimnames(gerb[[1]][[1]])[[2]]
nams.out <- gerb$nams.out
vars.new <- vars.old <- NULL
for (i in 1:length(vars)) {
if (is.element(vars[i], latent.vars)) {
tmp <- vars[i]
tmp1 <- nams.out[names(nams.out) == vars[i]]
} else if (is.element(vars[i], orig.vars)) {
tmp <- names(nams.out)[nams.out == vars[i]]
tmp1 <- rep(vars[i], length(tmp))
} else {
message(paste0("Variable name ", vars[i], " is invalid and will not be used for plotting."), "\n", appendLF = FALSE)
tmp <- tmp1 <- NULL
}
if (length(vars.new) == 0 & length(tmp) > 0) {
vars.new <- tmp
vars.old <- tmp1
} else if (length(tmp) > 0) {
vars.new[(length(vars.new) + 1):(length(vars.new) + length(tmp))] <- tmp
vars.old[(length(vars.old) + 1):(length(vars.old) + length(tmp1))] <- tmp1
}
}
if (length(vars.new) == 0) {
stop("No valid variable names were given.")
}
vars <- vars.new
if (is.logical(imps)) {
if (imps) {
imps <- as.numeric(dimnames(gerb.tmp)[[3]])
} else {
imps <- 1
}
}
if(length(imps) == 0) {
imps <- 1
}
if(length(obs.lwd) == 0) {
if (!sep & length(file) != 0) {
obs.lwd <- 1
} else {
obs.lwd <- 2
}
}
if(length(imp.lwd) == 0) {
if (!sep & length(file) != 0) {
imp.lwd <- 1
} else {
imp.lwd <- 2
}
}
if (length(imps) == 1) {
ltys <- imp.lty
lwds <- imp.lwd
} else {
if(length(imp.col) == 1) {
imp.col <- rep(imp.col, length(imps))
imp.col <- imp.col[1:length(imps)]
}
if(length(imp.lwd) == 1) {
imp.lwd <- rep(imp.lwd, length(imps))
imp.lwd <- imp.lwd[1:length(imps)]
}
if(length(imp.lty) == 1) {
imp.lty <- c(imp.lty, 3:(length(imps) + 1))
} else {
imp.lty <- rep(imp.lty, length(imps))
imp.lty <- imp.lty[1:length(imps)]
}
ltys <- imp.lty
lwds <- imp.lwd
}
for (i in 1:length(vars)) {
n.obs <- sum(gerb[[2]][,vars.old[i]] == 0)
n.mis <- sum(gerb[[2]][,vars.old[i]] != 0)
cont <- TRUE
if (n.mis == 0) {
cont <- FALSE
#warning(paste0("Will not produce univariate plot for ", var, ": No imputed values for the variable."))
message(paste0("Will not produce trace plots for ", vars[i], ": No imputed values for the variable."), "\n", appendLF = FALSE)
} else {
cont <- TRUE
}
if (cont) {
if (plot.obs) {
z <- gerb.tmp.obs[vars[i]]
} else {
z <- NULL
}
ylim <- c(min(c(gerb.tmp[vars[i], , imps], z)),max(c(gerb.tmp[vars[i], , imps], z)))
y <- gerb.tmp[vars[i], , imps[1]]
x <- as.numeric(dimnames(gerb.tmp)[[2]])
if (sep & length(file) > 0) {
if (!supp.suff) {
suff <- paste0("_", vars[i], "_time_series")
} else {
suff <- ""
}
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, suff, ".pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, suff, ".png", sep = ""),
width = width, height = height, units = "in", res = 500)
}
}
if (length(main) == 0) {
main1 <- vars[i]
} else {
main1 <- main
}
if (length(ylab) == 0) {
ylab1 <- vars[i]
} else {
ylab1 <- ylab
}
if (length(xlab) == 0) {
xlab1 <- "Iteration"
} else {
xlab1 <- xlab
}
plot(x, y, type = "l", lty = ltys[1], col = imp.col[1], main = main1, ylab = ylab1,
xlab = xlab1, ylim = ylim, lwd = lwds[1], panel.first = grid(), ...)
if (length(imps) > 1) {
for(j in 2:length(imps)) {
y <- gerb.tmp[vars[i], , imps[j]]
lines(x, y, lty = ltys[j], col = imp.col[j], lwd = lwds[j], ...)
}
}
leg <- paste0("Imputation ",imps)
cols <- imp.col
if (plot.obs) {
#obs.col <- rgb(0.3,0.2,0.4,0.6)
abline(h = z, col = obs.col, lty = obs.lty, lwd = obs.lwd)
leg <- c(leg, "Observed")
cols <- c(cols, obs.col)
ltys <- c(ltys, obs.lty)
lwds <- c(lwds, obs.lwd)
}
if (length(bty) == 0) {
bty <- "n"
}
if (plot.legend) {
if (length(legend.text) > 0) {
leg[1:length(legend.text)] <- legend.text
}
legend(legend.spot, legend = leg, bty = bty, col = cols, lty = ltys, lwd = lwds, ...)
}
if (sep & length(file) > 0) {
grDevices::dev.off()
}
}
}
if (!sep & length(file) > 0) {
grDevices::dev.off()
}
}
plot_gerbil.uni <- function (gerb, file = NULL, vars = NULL, imp = 1, obs.col = "blue4", imp.col = "brown2", sep = FALSE, obs.lty = 1, imp.lty = 2,
use.log = NULL, legend.text = NULL, legend.spot = "bottomleft", height = NULL, width = NULL, mfrow = c(3, 2), obs.lwd = NULL, imp.lwd = NULL,
main = NULL, xlab = NULL, ylab = NULL, bty = NULL, cex = NULL, ...) {
oldpar <- par(no.readonly = TRUE)
on.exit(par(oldpar))
plot.legend <- TRUE
if (is.logical(legend.text)) {
if(!legend) {
plot.legend <- FALSE
} else {
legend.text <- NULL
}
} else if (length(legend.text) == 1) {
if (legend.text == "n") {
plot.legend <- FALSE
} else if (legend.text == "y") {
legend.text <- NULL
}
}
if(length(obs.lwd) == 0) {
if (!sep & length(file) > 0) {
obs.lwd <- 1
} else {
obs.lwd <- 2
}
}
if(length(imp.lwd) == 0) {
if (!sep & length(file) > 0) {
imp.lwd <- 1
} else {
imp.lwd <- 2
}
}
if (is.logical(imp)) {
if (imp) {
imp <- as.numeric(dimnames(gerb)[[3]])
} else {
imp <- 1
}
}
if(length(imp) == 0) {
imp <- 1
}
#imp <- imp[1]
if(length(imp.col) == 1) {
imp.col <- rep(imp.col, length(imp))
imp.col <- imp.col[1:length(imp)]
}
if(length(imp.lwd) == 1) {
imp.lwd <- rep(imp.lwd, length(imp))
imp.lwd <- imp.lwd[1:length(imp)]
}
if(length(imp.lty) == 1) {
imp.lty <- c(imp.lty, 3:(length(imp) + 1))
} else {
imp.lty <- rep(imp.lty, length(imp))
imp.lty <- imp.lty[1:length(imp)]
}
col <- c(obs.col, imp.col)
lty <- c(obs.lty, imp.lty)
lwd <- c(obs.lwd, imp.lwd)
if (length(vars) == 0) {
#vars <- rownames(gerb[[3]])
#vars <- rownames(gerb$summary)[gerb$summary[, 4] > 0]
vars <- colnames(gerb$missing)[colSums(gerb$missing == 1 | gerb$missing == 2) > 0]
} else if (is.numeric(vars)) {
vars <- rownames(gerb[[3]])[vars]
}
if (length(file) == 0 | length(vars) == 1) {
sep <- TRUE
supp.suff <- TRUE
} else {
supp.suff <- FALSE
}
if (sep) {
mfrow <- c(1, 1)
}
if(length(height) == 0) {
if (!sep) {
height <- 11
} else {
height <- 11
}
}
if(length(width) == 0) {
if (!sep) {
width <- 8.5
} else {
width <- 8.5
}
}
if (length(file) == 0) {
graphics::par(mfrow = c(1, 1), 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 = mfrow, ask = FALSE)
}
}
if (is.logical(use.log)) {
if(use.log) {
use.log <- vars
}
}
for(i in 1:length(vars)) {
#type <- as.character(gerb[[3]][c(i), "Variable.Type"])
#var <- rownames(gerb[[3]])[c(i)]
type <- as.character(gerb[[3]][vars[i], "Variable.Type"])
var <- vars[i]
if (length(imp) > 1) {
dat <- matrix(NA, NROW(gerb[[1]][[imp[1]]]), length(imp))
for(j in 1:length(imp)) {
dat[, j] <- gerb[[1]][[imp[j]]][, var]
if (length(gerb$ineligibles) > 0) {
dat[gerb$ineligibles[, var], j] <- NA
}
}
} else {
dat <- gerb[[1]][[imp]][, var]
}
obs <- gerb[[2]][, var]
if(is.element(var, use.log)) {
if(sum(dat <= 0, na.rm = TRUE) > 0) {
dat[dat <= 0] <- NA
#warning(paste0("Non-positive values detected for ", var, ": Cases will be removed prior to log transformation."))
message(paste0("Non-positive values detected for ", var, ": Cases will be removed prior to log transformation."), "\n", appendLF = FALSE)
}
dat <- log(dat)
var.nam <- paste0("log(",var,")")
} else {
var.nam <- var
}
cens <- gerb[[2]][, var] == 2
if(sum(cens, na.rm = TRUE) > 0) {
#obs <- obs[!cens]
#dat <- dat[!cens]
if(length(imp) > 1) {
dat[cens, ] <- NA
} else {
dat[cens] <- NA
}
#warning(paste0("Censored values detected when plotting ", var, ": Censored cases will be removed."))
message(paste0("Censored values detected when plotting ", var, ": Censored cases will be removed."), "\n", appendLF = FALSE)
}
if (sum(obs == 1, na.rm = TRUE) == 0) {
cont <- FALSE
#warning(paste0("Will not produce univariate plot for ", var, ": No imputed values for the variable."))
message(paste0("Will not produce univariate plot for ", var, ": No imputed values for the variable."), "\n", appendLF = FALSE)
} else {
cont <- TRUE
}
#leg <- c("Observed", "Imputed")
if (length(imp) == 1) {
leg <- c("Observed", "Imputed")
} else {
leg <- c("Observed", paste0("Imputed: ", imp))
}
ns <- table(obs[!is.na(as.matrix(dat)[, 1])])
#leg <- paste0(leg, ": n = ", ns[c("0", "1")])
type <- as.character(gerb[[3]][var, "Variable.Type"])
if ((type == "binary" | type == "categorical" | type == "ordinal") & cont) {
if (length(file) == 0) {
par(mfrow = c(1, 1))
}
tab <- table(as.matrix(dat)[, 1], obs)
colnames(tab)[colnames(tab) == "1"] <- as.character(imp[1])
if (NCOL(dat) > 1) {
tab.tmp <- tab
tab <- matrix(NA, NROW(tab), 1 + NCOL(dat))
rownames(tab) <- rownames(tab.tmp)
colnames(tab) <- c("0", as.character(imp))
tab[, colnames(tab.tmp)] <- tab.tmp
for (j in 2:length(imp)) {
tab[, as.character(imp[j])] <- table(dat[, j], obs)[, "1"]
}
}
ests <- 100 * t(tab)/colSums(tab)
rownames(ests) <- leg
#ests <- 100 * rbind(
# Observed = tab[, "0"]/sum(tab[, "0"]),
# Imputed = tab[, "1"]/sum(tab[, "1"])
#)
#if (ests[1] <= ests[length(ests)]) {
# loc <- "topleft"
#} else {
# loc <- "topright"
#}
loc <- legend.spot
if (sep & length(file) > 0) {
if (!supp.suff) {
suff <- paste0("_", var)
} else {
suff <- ""
}
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, suff, ".pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, suff, ".png", sep = ""),
width = width, height = height, units = "in", res = 500)
}
}
if (length(main) == 0) {
main1 <- var.nam
} else {
main1 <- main
}
if (length(ylab) == 0) {
ylab1 <- "Frequency (%)"
} else {
ylab1 <- ylab
}
if (length(xlab) == 0) {
xlab1 <- ""
} else {
xlab1 <- xlab
}
if (length(bty) == 0) {
bty1 <- "n"
} else {
bty1 <- bty
}
if (plot.legend) {
if (length(legend.text) > 0) {
leg[1:length(legend.text)] <- legend.text
}
barplot(ests, beside = T, legend.text = T, col = col, args.legend = list(legend = leg, bty = bty1, x = loc), main = main1, ylab = ylab1, xlab = xlab1, ...)
} else {
barplot(ests, beside = T, col = col, main = main1, ylab = ylab1, xlab = xlab1, ...)
}
#grid()
abline(h = 0)
if (sep & length(file) > 0) {
grDevices::dev.off()
}
} else if (type == "semicont" & cont) {
mass <- gerb$mass.final[var]
if (length(file) == 0) {
#Edited by Max from (c(1, 2))
par(mfrow = c(2, 1))
}
#pos <- as.numeric(dat != mass)
#tab <- table(pos, obs)
#ests <- 100 * rbind(
# Observed = tab[, "0"]/sum(tab[, "0"]),
# Imputed = tab[, "1"]/sum(tab[, "1"])
#)
#colnames(ests) = paste0(c("Equals ","Does not equal "),mass)
#if (ests[1] <= ests[length(ests)]) {
# loc <- "topleft"
#} else {
# loc <- "topright"
#}
pos <- as.numeric(as.matrix(dat)[, 1] != mass)
tab <- table(pos, obs)
colnames(tab)[colnames(tab) == "1"] <- as.character(imp[1])
if (NCOL(dat) > 1) {
tab.tmp <- tab
tab <- matrix(NA, NROW(tab), 1 + NCOL(dat))
rownames(tab) <- rownames(tab.tmp)
colnames(tab) <- c("0", as.character(imp))
tab[, colnames(tab.tmp)] <- tab.tmp
for (j in 2:length(imp)) {
pos <- as.numeric(dat[, j] != mass)
tab[, as.character(imp[j])] <- table(pos, obs)[, "1"]
}
}
ests <- 100 * t(tab)/colSums(tab)
rownames(ests) <- leg
colnames(ests) <- paste0(c("Equals ","Does not equal "),mass)
loc <- legend.spot
if (sep & length(file) > 0) {
if (!supp.suff) {
suff <- paste0("_", var)
} else {
suff <- ""
}
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, suff, "_binary.pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, suff, "_binary.png", sep = ""),
width = width, height = height, units = "in", res = 500)
}
}
if (length(main) == 0) {
main1 <- paste0(var.nam, ": Binary")
} else {
main1 <- main
}
if (length(ylab) == 0) {
ylab1 <- "Frequency (%)"
} else {
ylab1 <- ylab
}
if (length(xlab) == 0) {
xlab1 <- ""
} else {
xlab1 <- xlab
}
if (length(bty) == 0) {
bty1 <- "n"
} else {
bty1 <- bty
}
if (length(cex) == 0) {
#cex1 <- .8
cex1 <- 1
} else {
cex1 <- cex
}
if (plot.legend) {
if (length(legend.text) > 0) {
leg[1:length(legend.text)] <- legend.text
}
barplot(ests, beside=T, legend.text = T, col = col,
args.legend = list(legend = leg, bty = bty1, x = loc, cex = cex1), main = main1, ylab = ylab1, xlab = xlab1, ...)
} else {
barplot(ests, beside=T, col = col, main = main1, ylab = ylab1, xlab = xlab1, ...)
}
#grid()
abline(h=0)
if (sep & length(file) > 0) {
grDevices::dev.off()
}
dat[dat == mass] <- NA
dat <- as.matrix(dat)
d1 <- dat[obs == 0, 1]
d1 <- d1[!is.na(d1)]
d1 <- density(d1)
if(is.element(var, use.log)) {
d1x <- exp(d1$x)
log.plot <- "x"
} else {
d1x <- d1$x
log.plot <- ""
}
d1y <- d1$y
xlim <- c(min(d1x), max(d1x))
ylim <- c(min(d1y), max(d1y))
d2y <- d2x <- list()
for(j in 1:length(imp)) {
d2 <- dat[obs == 1, j]
d2 <- d2[!is.na(d2)]
d2 <- density(d2)
if(is.element(var, use.log)) {
d2x[[j]] <- exp(d2$x)
} else {
d2x[[j]] <- d2$x
}
d2y[[j]] <- d2$y
xlim.tmp <- c(min(d2x[[j]]), max(d2x[[j]]))
ylim.tmp <- c(min(d2y[[j]]), max(d2y[[j]]))
xlim <- c(min(xlim[1], xlim.tmp[1]), max(xlim[2], xlim.tmp[2]))
ylim <- c(min(ylim[1], ylim.tmp[1]), max(ylim[2], ylim.tmp[2]))
}
#leg <- c("Observed","Imputed")
#ns <- c(length(d1),length(d2))
#leg <- paste0(leg,": n = ",ns)
if (sep & length(file) > 0) {
if (!supp.suff) {
suff <- paste0("_", var)
} else {
suff <- ""
}
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, suff, "_continuous.pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, suff, "_continuous.png", sep = ""),
width = width, height = height, units = "in", res = 500)
}
}
if (length(main) == 0) {
main1 <- paste0(var.nam, ": Continuous")
} else {
main1 <- main
}
if (length(ylab) == 0) {
ylab1 <- "Density"
} else {
ylab1 <- ylab
}
if (length(xlab) == 0) {
xlab1 <- ""
} else {
xlab1 <- xlab
}
plot(d1x, d1y, xlim = xlim, ylim = ylim, lty = lty[1], col = col[1], type = "l", xlab = xlab1, ylab = ylab1,
main = main1, lwd = lwd[1], log = log.plot, panel.first = grid(), ...)
for(j in 1:length(imp)) {
lines(d2x[[j]], d2y[[j]], xlim = xlim, ylim = ylim, lty = lty[j + 1], col = col[j + 1], type = "l", lwd = lwd[j + 1])
}
#Edited by Max below, to include the cex option
if (plot.legend) {
if (length(legend.text) > 0) {
leg[1:length(legend.text)] <- legend.text
}
legend(legend.spot, legend = leg, bty = bty1, col = col, lty = lty, lwd = lwd, cex = cex1, ...)
}
if (sep & length(file) > 0) {
grDevices::dev.off()
}
} else if (cont) {
if (length(file) == 0) {
par(mfrow = c(1, 1))
}
dat <- as.matrix(dat)
d1 <- dat[obs == 0, 1]
d1 <- d1[!is.na(d1)]
d1 <- density(d1)
if(is.element(var, use.log)) {
d1x <- exp(d1$x)
log.plot <- "x"
} else {
d1x <- d1$x
log.plot <- ""
}
d1y <- d1$y
xlim <- c(min(d1x), max(d1x))
ylim <- c(min(d1y), max(d1y))
d2y <- d2x <- list()
for(j in 1:length(imp)) {
d2 <- dat[obs == 1, j]
d2 <- d2[!is.na(d2)]
d2 <- density(d2)
if(is.element(var, use.log)) {
d2x[[j]] <- exp(d2$x)
} else {
d2x[[j]] <- d2$x
}
d2y[[j]] <- d2$y
xlim.tmp <- c(min(d2x[[j]]), max(d2x[[j]]))
ylim.tmp <- c(min(d2y[[j]]), max(d2y[[j]]))
xlim <- c(min(xlim[1], xlim.tmp[1]), max(xlim[2], xlim.tmp[2]))
ylim <- c(min(ylim[1], ylim.tmp[1]), max(ylim[2], ylim.tmp[2]))
}
#lty <- c(imp.lty, obs.lty)
#lty <- c(obs.lty, ltys)
#col <- c(1, 2)
if (sep & length(file) > 0) {
if (!supp.suff) {
suff <- paste0("_", var)
} else {
suff <- ""
}
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, suff, ".pdf", sep = ""), width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, suff, ".png", sep = ""),
width = width, height = height, units = "in", res = 500)
}
}
if (length(main) == 0) {
main1 <- var.nam
} else {
main1 <- main
}
if (length(ylab) == 0) {
ylab1 <- "Density"
} else {
ylab1 <- ylab
}
if (length(bty) == 0) {
bty1 <- "n"
} else {
bty1 <- bty
}
if (length(cex) == 0) {
cex1 <- .8
} else {
cex1 <- cex
}
if (length(xlab) == 0) {
xlab1 <- ""
} else {
xlab1 <- xlab
}
plot(d1x, d1y, xlim = xlim, ylim = ylim, lty = lty[1], col = col[1], type = "l", xlab = xlab1,
ylab = ylab1, main = main1, lwd = lwd[1], log = log.plot, panel.first = grid(), ...)
for(j in 1:length(imp)) {
lines(d2x[[j]], d2y[[j]], xlim = xlim, ylim = ylim, lty = lty[j + 1], col = col[j + 1], type = "l", lwd = lwd[j + 1])
}
if (plot.legend) {
if (length(legend.text) > 0) {
leg[1:length(legend.text)] <- legend.text
}
legend(legend.spot, legend = leg, bty = bty1, col = col, lty = lty, lwd = lwd, ...)
}
if (sep & length(file) > 0) {
grDevices::dev.off()
}
}
}
if (!sep & length(file) > 0) {
grDevices::dev.off()
}
}
plot_gerbil.mult <- function(gerb, file = NULL, vars = NULL, imp = 1, obs.col = "blue4", imp.col = "brown2", sep = FALSE, legend.text = NULL,
legend.spot = "bottomleft", height = NULL, width = NULL, mfrow = c(3, 2), pch = c(1, 3), use.log = NULL,
lwd = NULL, main = NULL, xlab = NULL, ylab = NULL, bty = NULL, cex = NULL, lty = NULL, xaxt = "s", yaxt = "s",
partial = "imputed", ...) {
oldpar <- par(no.readonly = TRUE)
on.exit(par(oldpar))
plot.legend <- TRUE
if (is.logical(legend.text)) {
if(!legend) {
plot.legend <- FALSE
} else {
legend.text <- NULL
}
} else if (length(legend.text) == 1) {
if (legend.text == "n") {
plot.legend <- FALSE
} else if (legend.text == "y") {
legend.text <- NULL
}
}
if (is.logical(imp)) {
if (imp) {
imp <- as.numeric(dimnames(gerb)[[3]])
} else {
imp <- 1
}
}
if(length(imp) == 0) {
imp <- 1
}
if(length(imp.col) == 1) {
imp.col <- rep(imp.col, length(imp))
imp.col <- imp.col[1:length(imp)]
}
col <- c(obs.col, imp.col)
obs.pch <- pch[1]
if (length(pch) == 1) {
imp.pch <- pch[1]
} else {
imp.pch <- pch[-1]
}
if(length(imp.pch) == 1) {
imp.pch <- rep(imp.pch, length(imp))
imp.pch <- imp.pch[1:length(imp)]
}
pch <- c(obs.pch, imp.pch)
obs.lwd <- lwd[1]
if (length(lwd) == 1) {
imp.lwd <- lwd[1]
} else {
imp.lwd <- lwd[-1]
}
if(length(imp.lwd) == 1) {
imp.lwd <- rep(imp.lwd, length(imp))
imp.lwd <- imp.lwd[1:length(imp)]
}
lwd <- c(obs.lwd, imp.lwd)
obs.lty <- lty[1]
if (length(lty) == 1) {
imp.lty <- lty[1]
} else {
imp.lty <- lty[-1]
}
if(length(imp.lty) == 1) {
imp.lty <- rep(imp.lty, length(imp))
imp.lty <- imp.lty[1:length(imp)]
}
lty <- c(obs.lty, imp.lty)
#imp <- imp[1]
if (is.logical(vars)) {
if (vars) {
vars <- colnames(gerb$missing)[colSums(gerb$missing == 1 | gerb$missing == 2) > 0]
}
else {
vars <- colnames(gerb$missing)
}
}
else if (length(vars) == 0) {
vars <- colnames(gerb$missing)
}
if (length(file) == 0 | length(vars) == 2) {
sep <- TRUE
supp.suff <- TRUE
} else {
supp.suff <- FALSE
}
if (sep) {
mfrow <- c(1, 1)
}
if(length(vars) <= 1) {
stop("Need at least 2 variables for bivariate plotting.")
}
if(length(height) == 0) {
if (!sep) {
height <- 11
} else {
height <- 6.5
}
}
if(length(width) == 0) {
if (!sep) {
width <- 8.5
} else {
width <- 6.5
}
}
if (is.logical(use.log)) {
if(use.log) {
use.log <- vars
}
}
if (length(file) == 0) {
#graphics::par(mfrow = c(1, 1), 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 = mfrow, ask = FALSE)
}
}
if (length(imp) == 1) {
leg <- c("Observed", "Imputed")
leg <- c("Obs.", "Imp.")
} else {
leg <- c("Observed", paste0("Imputed ", imp))
}
#ns <- table(obs[!is.na(as.matrix(dat)[, 1])])
if (length(legend.text) > 0) {
leg[1:length(legend.text)] <- legend.text
}
plot.num <- 0
neg.vars <- NA
for(j in 2:length(vars)) {
for(i in 1:(j-1)) {
types <- gerb[[3]][c(vars[i], vars[j]), "Variable.Type"]
#vars1 <- rownames(gerb[[3]])[c(i, j)]
vars1 <- c(vars[i], vars[j])
is.box <- types == "categorical" | types == "ordinal" | types == "binary"
if(!is.box[1] & is.box[2]) {
var1 <- vars1[2]
var2 <- vars1[1]
} else if (is.box[1] & is.box[2]) {
dat1 <- gerb[[1]][[imp[1]]][, vars1[1]]
dat2 <- gerb[[1]][[imp[1]]][, vars1[2]]
if (length(gerb$ineligibles) > 0) {
dat1[gerb$ineligibles[, vars1[1]]] <- NA
dat2[gerb$ineligibles[, vars1[2]]] <- NA
}
#n1 <- length(table(gerb[[1]][[imp]][, vars1[1]]))
#n2 <- length(table(gerb[[1]][[imp]][, vars1[2]]))
n1 <- length(table(dat1))
n2 <- length(table(dat2))
if(n2 > n1) {
var2 <- vars1[2]
var1 <- vars1[1]
} else {
var2 <- vars1[1]
var1 <- vars1[2]
}
} else {
var1 <- vars1[1]
var2 <- vars1[2]
}
if (partial == "imputed") {
miss <- gerb[[2]][, var1] == 1 | gerb[[2]][, var2] == 1
} else {
miss <- gerb[[2]][, var1] == 1 & gerb[[2]][, var2] == 1
}
if (length(imp) > 1) {
dat1 <- dat2 <- matrix(NA, NROW(gerb[[1]][[imp[1]]]), length(imp))
for(k in 1:length(imp)) {
c.d1 <- class(gerb[[1]][[imp[k]]][, var1])[1]
c.d2 <- class(gerb[[1]][[imp[k]]][, var2])[1]
if (c.d1 == "character" | c.d1 == "factor" | c.d1 == "ordered") {
d1.tmp <- as.character(gerb[[1]][[imp[k]]][, var1])
} else {
d1.tmp <- gerb[[1]][[imp[k]]][, var1]
}
if (c.d2 == "character" | c.d2 == "factor" | c.d2 == "ordered") {
d2.tmp <- as.character(gerb[[1]][[imp[k]]][, var2])
} else {
d2.tmp <- gerb[[1]][[imp[k]]][, var2]
}
dat1[, k] <- d1.tmp
dat2[, k] <- d2.tmp
if (length(gerb$ineligibles) > 0) {
dat1[gerb$ineligibles[, var1], k] <- NA
dat2[gerb$ineligibles[, var2], k] <- NA
}
}
} else {
c.d1 <- class(gerb[[1]][[imp]][, var1])[1]
c.d2 <- class(gerb[[1]][[imp]][, var2])[1]
if (c.d1 == "character" | c.d1 == "factor" | c.d1 == "ordered") {
d1.tmp <- as.character(gerb[[1]][[imp]][, var1])
} else {
d1.tmp <- gerb[[1]][[imp]][, var1]
}
if (c.d2 == "character" | c.d2 == "factor" | c.d2 == "ordered") {
d2.tmp <- as.character(gerb[[1]][[imp]][, var2])
} else {
d2.tmp <- gerb[[1]][[imp]][, var2]
}
dat1 <- d1.tmp
dat2 <- d2.tmp
if (length(gerb$ineligibles) > 0) {
dat1[gerb$ineligibles[, var1]] <- NA
dat2[gerb$ineligibles[, var2]] <- NA
}
}
cens <- gerb[[2]][, var1] == 2 | gerb[[2]][, var2] == 2
if(sum(cens, na.rm = TRUE) > 0) {
#miss <- miss[!cens]
#dat1 <- dat1[!cens]
#dat2 <- dat2[!cens]
if(length(imp) > 1) {
dat1[cens, ] <- NA
dat2[cens, ] <- NA
} else {
dat1[cens] <- NA
dat2[cens] <- NA
}
#warning(paste0("Censored values detected when plotting ", var1, " and ", var2, ": Censored cases will be removed."))
message(paste0("Censored values detected when plotting ", var1, " and ", var2, ": Censored cases will be removed."), "\n", appendLF = FALSE)
}
if (sum(miss, na.rm = TRUE) == 0) {
cont <- FALSE
#warning(paste0("Will not produce bivariate plot between ", var1, " and ", var2, ": No imputed values for either variable."))
message(paste0("Will not produce bivariate plot between ", var1, " and ", var2, ": No imputed values for either variable."), "\n", appendLF = FALSE)
} else {
cont <- TRUE
}
if(is.element(var1, use.log)) {
if (sum(dat1 <= 0, na.rm = TRUE) > 0) {
dat1[dat1 <= 0] <- NA
if (!is.element(var1, neg.vars)) {
#warning(paste0("Non-positive values detected for ", var1, ": Cases will be removed prior to log transformation."))
message(paste0("Non-positive values detected for ", var1, ": Cases will be removed prior to log transformation."), "\n", appendLF = FALSE)
neg.vars[(length(neg.vars) + 1)] <- var1
}
}
dat1 <- log10(dat1)
var.nam1 <- paste0("log(",var1,")")
} else {
var.nam1 <- var1
}
if(is.element(var2, use.log)) {
if (sum(dat2 <= 0, na.rm = TRUE) > 0) {
dat2[dat2 <= 0] <- NA
if (!is.element(var2, neg.vars)) {
#warning(paste0("Non-positive values detected for ", var2, ": Cases will be removed prior to log transformation."))
message(paste0("Non-positive values detected for ", var2, ": Cases will be removed prior to log transformation."), "\n", appendLF = FALSE)
neg.vars[(length(neg.vars) + 1)] <- var2
}
}
dat2 <- log10(dat2)
var.nam2 <- paste0("log(",var2,")")
} else {
var.nam2 <- var2
}
if (sum(is.box) == 2) {
#warning(paste0("Will not produce bivariate plot between ", var1, " and ", var2, ": Neither variable is continuous."))
#message(paste0("Will not produce bivariate plot between ", var1, " and ", var2, ": Neither variable is continuous."), "\n", appendLF = FALSE)
#cont <- FALSE
}
#leg <- c("Observed", "Imputed")
#leg <- paste0(leg, ": n = ", ns[c("0", "1")])
if (!supp.suff) {
suff <- paste0("_", var2, "_vs_", var1)
} else {
suff <- ""
}
if (sep & length(file) > 0 & cont) {
if (file.type == "pdf") {
grDevices::pdf(file = paste(file, suff, ".pdf", sep = ""),
width = width, height = height)
} else if (file.type == "png") {
grDevices::png(file = paste(file, suff, ".png", sep = ""),
width = width, height = height, units = "in", res = 500)
}
}
if (sum(is.box) == 0) {
num.plots <- length(imp)
} else {
num.plots <- 1
}
plot.num <- plot.num + cont * num.plots
if(sum(is.box) == 0 & cont) {
ylim <- c(min(dat2, na.rm = TRUE), max(dat2, na.rm = TRUE))
xlim <- c(min(dat1, na.rm = TRUE), max(dat1, na.rm = TRUE))
#pch <- c("o","x")
#col <- c(2,1)
if (length(ylab) == 0) {
#ylab1 <- var.nam2
ylab1 <- var2
} else {
ylab1 <- ylab
}
if (length(bty) == 0) {
bty1 <- "n"
} else {
bty1 <- bty
}
if (length(xlab) == 0) {
#xlab1 <- var.nam1
xlab1 <- var1
} else {
xlab1 <- xlab
}
dat1 <- as.matrix(dat1)
dat2 <- as.matrix(dat2)
if(is.element(var1, use.log)) {
x1 <- exp(dat1[!miss, ])
x2 <- exp(dat1[miss, ])
log.plot <- "x"
} else {
x1 <- dat1[!miss, ]
x2 <- dat1[miss, ]
log.plot <- ""
}
if(is.element(var2, use.log)) {
y1 <- exp(dat2[!miss, ])
y2 <- exp(dat2[miss, ])
log.plot <- paste0(log.plot, "y")
} else {
y1 <- dat2[!miss, ]
y2 <- dat2[miss, ]
}
if (length(main) == 0) {
main1 <- paste0(var.nam2," vs. ", var.nam1)
} else {
main1 <- main
}
if (NCOL(dat1) > 1) {
message(paste0("Both variables ", var1, " and ", var2, " are continuous: Separate plots will be made for each element of 'imp'."), "\n", appendLF = FALSE)
}
for(k in 1:NCOL(dat1)) {
if (NCOL(dat1) > 1) {
main2 <- paste0(leg[k + 1], ": ", main1)
} else {
main2 <- main1
}
#plot(dat1[!miss], dat2[!miss], pch = pch[1], col = col[1], xlim = xlim, ylim = ylim, main = main2, xlab = xlab1, ylab = ylab1, ...)
plot(dat1[!miss, k], dat2[!miss, k], pch = pch[1], col = col[1], xlim = xlim, ylim = ylim, main = main2, xlab = xlab1, ylab = ylab1, xaxt = "n", yaxt = "n", ...)
ats1 <- par("xaxp")
ats2 <- par("yaxp")
ats1 <- seq(ats1[1], ats1[2], length.out = ats1[3] + 1)
if(is.element(var1, use.log)) {
labs1 <- 10^ats1
} else {
labs1 <- ats1
}
if (is.element(xaxt, c("s", "l", "t"))) {
axis(side = 1, at = ats1, labels = labs1)
}
ats2 <- seq(ats2[1], ats2[2], length.out = ats2[3] + 1)
if(is.element(var2, use.log)) {
labs2 <- 10^ats2
} else {
labs2 <- ats2
}
if (is.element(yaxt, c("s", "l", "t"))) {
axis(side = 2, at = ats2, labels = labs2)
}
abline(lm(dat2[!miss, k] ~ dat1[!miss, k]), lty = lty[1], col = col[1], lwd = lwd[1])
points(dat1[miss, k], dat2[miss, k], pch = pch[k + 1], col = col[k + 1])
abline(lm(dat2[miss, k] ~ dat1[miss, k]), lty = lty[k + 1], col = col[k + 1], lwd = lwd[k + 1])
if (plot.legend) {
legend(legend.spot, legend = leg[c(1, k + 1)], bty = bty1, col = col[c(1, k + 1)], pch = pch[c(1, k + 1)], ...)
}
}
} else if (sum(is.box) == 1 & cont) {
tmp.df <- data.frame(Y = as.matrix(dat2)[, 1], X = as.matrix(dat1)[, 1], R = miss)
if (NCOL(dat1) > 1) {
n.m <- sum(miss)
tmp.df1 <- matrix(NA, n.m * (NCOL(dat1) - 1), 3)
colnames(tmp.df1) <- colnames(tmp.df)
for (k in 2:NCOL(dat1)) {
tmp.df1[(((k - 1) - 1) * n.m + 1):((k - 1) * n.m), ] <- cbind(Y = dat2[miss, k], X = dat1[miss, k], R = k)
}
tmp.df[, "R"] <- as.numeric(tmp.df[, "R"])
tmp.df1 <- data.frame(tmp.df1)
tmp.df <- rbind(tmp.df, tmp.df1)
}
tmp.df <- tmp.df[!is.na(tmp.df$Y) & !is.na(tmp.df$X),]
tmp.df$Y <- as.numeric(as.character(tmp.df$Y))
tab <- table(tmp.df$X)
p <- length(tab)
p1 <- NCOL(dat1) + 1
sep1 <- 1
at <- 1:(p1 * p) + c(t(matrix(rep(sep1 * (0:(p - 1)), p1), p, p1)))
here <- matrix(at, p1, p)
here <- colMeans(here)
if (length(main) == 0) {
main1 <- paste0(var.nam2, " vs. ", var.nam1)
} else {
main1 <- main
}
if (length(ylab) == 0) {
#ylab1 <- var.nam2
ylab1 <- var2
} else {
ylab1 <- ylab
}
if (length(bty) == 0) {
bty1 <- "n"
} else {
bty1 <- bty
}
if (length(xlab) == 0) {
xlab1 <- var.nam1
} else {
xlab1 <- xlab
}
boxplot(Y ~ R + X, data = tmp.df, at = at, col = col, xaxt = "n", yaxt = "n", ylab = ylab1, xlab = xlab1, main = main1, ...)
ats2 <- par("yaxp")
ats2 <- seq(ats2[1], ats2[2], length.out = ats2[3] + 1)
if(is.element(var2, use.log)) {
labs2 <- 10^ats2
} else {
labs2 <- ats2
}
if (is.element(yaxt, c("s", "l", "t"))) {
axis(side = 2, at = ats2, labels = labs2)
}
if (is.element(xaxt, c("s", "l", "t"))) {
axis(side = 1, at = here, labels = names(tab), tick = TRUE)
}
if (plot.legend) {
legend(legend.spot, legend = leg, fill = col, bty = bty1, ...)
}
} else if (sum(is.box) == 2 & cont) {
#miss <- gerb[[2]][, var1] == 1 | gerb[[2]][, var2] == 1
#dat1 <- gerb[[1]][[imp]][, var1]
#dat2 <- gerb[[1]][[imp]][, var2]
tab <- table(as.matrix(dat1)[, 1], as.matrix(dat2)[, 1], miss)
dimnames(tab)[[3]] <- c("0","1")
if (NCOL(dat1) > 1) {
p1 <- NCOL(dat1)
tab.tmp <- array(NA, c(dim(tab)[1], dim(tab)[2], p1 + 1))
dimnames(tab.tmp) <- list(dimnames(tab)[[1]], dimnames(tab)[[2]], 0:p1)
tab.tmp[, , 1:2] <- tab
for(k in 2:p1) {
tab.tmp1 <- table(dat1[, k], dat2[, k], miss)
tab.tmp[, , k + 1] <- tab.tmp1[, , "TRUE"]
}
tab <- tab.tmp
}
sums.tab <- apply(tab, 3, sum)
for (k in 1:dim(tab)[3]) {
tab[, , k] <- tab[, , k]/sums.tab[k]
}
tab <- 100 * tab
nams1 <- dimnames(tab)[[1]]
nams2 <- dimnames(tab)[[2]]
nams3 <- dimnames(tab)[[3]]
if(sum(is.na(as.numeric(nams1))) == 0) {
nams1 <- as.numeric(nams1)
}
if(sum(is.na(as.numeric(nams2))) == 0) {
nams2 <- as.numeric(nams2)
}
x1.tmp <- factor(rep(nams1, length(nams2)))
x2.tmp <- factor(c(t(matrix(nams2, length(nams2), length(nams1)))))
x3.tmp <- matrix(nams3, length(nams3), length(x1.tmp))
x3.tmp <- c(t(x3.tmp))
x1.tmp <- rep(x1.tmp, length(nams3))
x2.tmp <- rep(x2.tmp, length(nams3))
tmp <- data.frame(x1 = x1.tmp, x2 = x2.tmp, x3 = x3.tmp)
#tmp <- data.frame(
# x1 = factor(rep(nams1, length(nams2))),
# x2 = factor(c(t(matrix(nams2, length(nams2), length(nams1)))))
#)
#tmp <- rbind(data.frame(tmp, x3 = nams3[1]), data.frame(tmp, x3 = nams3[2]))
tmp1 <- rep(NA,NROW(tmp))
for(ind in 1:NROW(tmp)) {
tmp1[ind] <- tab[tmp[ind, 1], tmp[ind, 2], tmp[ind, 3]]
}
tmp[, 3] <- as.character(tmp[, 3])
#tmp[tmp[, 3] == "FALSE", 3] <- "Observed"
#tmp[tmp[, 3] == "TRUE", 3] <- "Imputed"
tmp[, 3] <- leg[match(tmp[, 3], nams3)]
colnames(tmp) <- c(var1, var2, "Status")
tmp2 <- data.frame(tmp, Frequency = tmp1)
tmp2$Status <- factor(tmp2$Status, levels = leg)
# require(lattice)
#pch = c(1, 3)
#pch <- c("o","x")
#cols=c("blue4", "brown2")
form <- paste0(var2, "~", "Frequency", "|", var1)
form <- as.formula(form)
if (length(file) == 0) {
size.tmp <- dev.size()
dev.off()
dev.new(height = size.tmp[2], width = size.tmp[1])
}
if(!sep) {
plot.new()
num.tmp <- plot.num %% prod(mfrow)
if(num.tmp == 0) {
num.tmp <- prod(mfrow)
}
n.y <- ceiling(num.tmp/mfrow[2])
n.x <- num.tmp - mfrow[2]*(n.y - 1)
if (length(cex) == 0) {
cex2 <- 0.7
cex1 <- 0.8
} else {
cex2 <- cex
cex1 <- cex
}
} else {
if (length(cex) == 0) {
cex2 <- 1
cex1 <- 1
} else {
cex2 <- cex
cex1 <- cex
}
n.x <- n.y <- 1
}
if (length(main) == 0) {
main1 <- paste0(var.nam2," vs. ", var.nam1)
} else {
main1 <- main
}
if (length(bty) == 0) {
bty1 <- "n"
} else {
bty1 <- bty
}
if (length(xlab) == 0) {
xlab1 <- "Frequency (%)"
} else {
xlab1 <- xlab
}
if (length(ylab) == 0) {
ylab1 <- ""
} else {
ylab1 <- ylab
}
if (plot.legend) {
Status <- tmp2$Status
plot(
lattice::dotplot(form,
data = tmp2,
horizontal = TRUE,
#panel = panel.superpose,
group = Status, pch = pch, col = col,
scales = list(cex = cex2),
par.strip.text = list(cex = cex2),
key = list(space = "right",
#transparent = TRUE,
points = list(pch = pch, col = col),
text = list(leg, cex = cex2),
cex.title = cex2, cex = cex2),
main = list(main1, cex = cex1), xlab = list(xlab1, cex = cex2), ylab = list(ylab1, cex = cex2), cex = cex2),
split = c(n.x, n.y, mfrow[2], mfrow[1]), cex = cex2,
newpage = FALSE
)
} else {
Status <- tmp2$Status
plot(
lattice::dotplot(form,
data = tmp2,
horizontal = TRUE,
#panel = panel.superpose,
group = Status, pch = pch, col = col,
scales = list(cex = cex2),
par.strip.text = list(cex = cex2),
main = list(main1, cex = cex1), xlab = list(xlab1, cex = cex2), ylab = list(ylab1, cex = cex2), cex = cex2),
split = c(n.x, n.y, mfrow[2], mfrow[1]), cex = cex2,
newpage = FALSE
)
}
}
if (sep & length(file) > 0 & cont) {
grDevices::dev.off()
}
}
}
if (!sep & length(file) > 0) {
grDevices::dev.off()
}
}
#' @title Plotting for gerbil objects
#'
#' @description
#' Using a \code{gerbil} object as an input, this function gives
#' diagnostic plots for selected variables
#'
#' @details
#' Three types of plots may be produced:
#' 1) Univariate (produced by setting \code{type = 1}): Compares the marginal distribution of observed and imputed values of a given variable. Density plots are produced for continuous variables, and bar plots are given for binary, categorical, and ordinal variables. For semi-continuous variables, two plots are constructed: a) a bar plot for the binary portion of the variable and 2) a density plot for the continuous portion.
#' 2) Bivariate (produced by setting \code{type = 2}): Compares the bivariate distributions of observed and imputed values of two variables. Scatter plots are produced if both variables are continuous or semi-continuous, box plots are produced if one variable is continuous or semi-continuous and the other is not, and a lattice plot is produced if neither variable is continuous or semi-continuous. For bivariate plots, imputed observations are those that have one or more of the values of the pair missing within the original dataset.
#' 3) Trace lines (produced by setting \code{type = 3}): Plots a pre-specified parameter across iterations of MCMC in order to examine convergence for a given variable. Parameters that may be plotted include means (\code{trace.type = 1}) and variances (\code{trace.type = 2}).
#'
#' Multiple plots may be created, as determined by the variable names listed in the parameter \code{y}. For univariate and trace plots, one plot is created for
#' each variable listed in \code{y}. For bivariate plotting, one plot is created for each combination of two elements within the vector \code{y} (as such, \code{y} must have a length of at least two in this case).
#' For trace plotting, elements of \code{y} should correspond to column names in the dataset that has been expanded to include binary indicators for categorical and semi-continuous variables.
#' If multiple plots are to be created, it is recommended to specify a file for output using the parameter \code{file}, in which case separate
#' files will be created for each plot (if \code{sep = TRUE}) or all plots will be written to the same file (if \code{sep = FALSE}).
#'
#' The only required input is a parameter \code{x} which is a \code{gerbil} object.
#'
#' @param x A \code{gerbil} object containing the imputed data.
#' @param y A vector listing the column names of the imputed data for which plots should be created. See details. By default, \code{y} contains all columns of the data that required imputation.
#' @param type A scalar used to specify the type of plots that will be created. Options include univariate (marginal) plots (\code{type = 1}), bivariate plots (\code{type = 2}), and trace plots (\code{type = 3}). See details. Defaults to \code{type = 1}.
#' @param imp A scalar or vector indicating which of the multiply imputed datasets should be used for plotting. Defaults to \code{imp = 1}. Setting \code{imp = TRUE} will include all imputed datasets.
#' @param col The color used for plotting -- should be a vector of length equal to \code{imp + 1}. The first element references plotting of observed data, and remaining elements reference plotting of imputed data.
#' @param lty The line type used for plotting imputed values with trace lines or density plots -- should be a vector of length equal to \code{imp + 1}. The first element references plotting of observed data, and remaining elements reference plotting of imputed data.
#' @param lwd The line width used for density and trace line plotting -- should be a vector of length equal to \code{imp + 1}. The first element references plotting of observed data, and remaining elements reference plotting of imputed data.
#' @param pch A length-2 vector that indicates the plotting symbol to be used for imputed and observed values in scatter and lattice plots.
#' @param log A character vector that includes names of variables of which a log transformation is to be taken prior to plotting.
#' @param legend A character or expression vector to appear in the legend. If \code{FALSE} or \code{'n'}, no legend is created. Defaults to \code{c("Observed", "Imputed", ...)}.
#' @param legend.loc The location of the legend in the plots.
#' @param mfrow The layout of plots across a single page when there are to be multiple plots per page (as is the case when \code{file} is non-\code{NULL} and \code{sep = FALSE}).
#' @param trace.type The type of trace plot to be created (only valid when \code{type = 3}). See details. Defaults to \code{trace.type = 1}.
#' @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} or \code{.png} extension. If \code{NULL} (the default), no file is created.
#' @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 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 partial Indicates how partially imputed pairs are handled in bivariate plotting. If \code{'imputed'}, cases with at least one missing variable in a pair are considered imputed. Otherwise (\code{partial = 'observed'}), only cases with both variables in the pair missing are considered imputed.
#' @param ... Arguments to be passed to methods, such as \code{plot}.
#'
#' @return No returned value, but instead plots are generated in the workspace or written to a specified directory.
#'
#' @examples
#' \donttest{
#' #Load the India Human Development Survey-II dataset
#' data(ihd_mcar)
#'
#' # Create a gerbil object
#'
#' imps.gerbil <- gerbil(ihd_mcar, m = 1, ords = "education_level", semi = "farm_labour_days",
#' bincat = "job_field")
#'
#' # Univariate plotting of all variables to a file
#' plot(imps.gerbil, type = 1, file = file.path(tempdir(), "gerbil_univariate.pdf"))
#'
#' # Bivariate plotting of all variables to a file
#' plot(imps.gerbil, type = 2, file = file.path(tempdir(), "gerbil_bivariate.pdf"))
#'
#' # Trace plotting of all variables to a file
#' plot(imps.gerbil, type = 3, file = file.path(tempdir(), "gerbil_ts.pdf"))
#'
#' # Univariate plotting of one variable (not to a file)
#' plot(imps.gerbil, type = 1, y = "job_field")
#'
#' # Bivariate plotting of one pair of variables (not to a file)
#' plot(imps.gerbil, type = 2, y = c("job_field", "income"))
#'
#' # Bivariate plotting of one pair of variables (not to a file) with income logged
#' plot(imps.gerbil, type = 2, y = c("job_field", "income"), log = "income")
#' }
#'
#' @export
#'
#'
#' @importFrom stats as.formula density
#' @importFrom graphics abline axis barplot boxplot grid legend lines par plot plot.new points
#' @importFrom grDevices dev.new dev.off dev.size
#'
plot.gerbil <- function(x, y = NULL, type = "Univariate", imp = 1, col = NULL, lty = NULL, lwd = NULL, pch = NULL,
log = NULL, legend = NULL, legend.loc = "topright", mfrow = c(3, 2), trace.type = "Mean",
file = NULL, sep = FALSE, height = NULL, width = NULL, partial = "imputed", ...)
{
if(length(file) == 0) {
#if ((tolower(type) == "univariate" | type == 1) & length(y) != 1) {
# warning("y should have length 1 for univariate plotting.")
#} else
if ((tolower(type) == "bivariate" | tolower(type) == "multivariate" | type == 2) & length(y) < 2) {
warning("y should have length of at least 2 for bivariate plotting.")
}
#else if ((tolower(type) == "ts" | tolower(type) == "time.series" | type == 3) & length(y) != 1) {
# warning("y should have length 1 for trace plotting.")
#}
}
if (length(pch) == 0) {
pch <- c(1, 3)
} else if (length(pch) == 1) {
pch <- c(pch, pch)
}
if (length(col) == 0) {
col <- c("blue4", "brown2")
} else if (length(col) == 1) {
col <- rep(col, 2)
}
if (length(lty) == 0) {
lty <- c(2, 1)
} else if (length(lty) == 1) {
lty <- rep(lty, 2)
}
obs.col <- col[1]
imp.col <- col[-1]
obs.lty <- lty[1]
imp.lty <- lty[-1]
if (length(lwd) <= 1) {
imp.lwd <- obs.lwd <- lwd
} else {
obs.lwd <- lwd[1]
imp.lwd <- lwd[-1]
}
plot_gerbil(gerb = x, vars = y, type = type,
imp = imp, trace.type = trace.type, log = log,
obs.col = obs.col, imp.col = imp.col, legend.text = legend,
legend.spot = legend.loc, obs.lty = obs.lty, imp.lty = imp.lty,
obs.lwd = obs.lwd, imp.lwd = imp.lwd, pch = pch,
file = file, sep = sep, mfrow = mfrow, height = height,
width = width, partial = partial, ...)
}
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.