On Thu, 18 Apr 2002, Mark Harris wrote:
> Hello,
>
> Can anyone help me out with this problem?
> After performing logistic regressions and testing the significance with
> likelihood ratios, I have plotted the results using "termplot". I
am
> wondering, how to get the names of my variables to appear on the x-axis
> rather than ascending numbers?
> I have used:
>
The following revised version of termplot() does this (not extensively
checked, but it works on a few examples). It's a two-line change but as
R1.5.0 is already frozen this won't get added until the next version.
-thomas
"termplot" <-
function (model, data = model.frame(model), partial.resid = FALSE,
rug = FALSE, terms = NULL, se = FALSE, xlabs = NULL, ylabs = NULL,
main = NULL, col.term = 2, lwd.term = 1.5, col.se = "orange",
lty.se = 2, lwd.se = 1, col.res = "gray", cex = 1, pch =
par("pch"),
ask = interactive() && nb.fig < n.tms && .Device !=
"postscript",
...)
{
terms <- if (is.null(terms))
predict(model, type = "terms", se = se)
else predict(model, type = "terms", se = se, terms = terms)
n.tms <- ncol(tms <- as.matrix(if (se)
terms$fit
else terms))
mf <- model.frame(model)
nmt <- colnames(tms)
cn <- parse(text = nmt)
if (is.null(ylabs))
ylabs <- paste("Partial for", nmt)
if (is.null(main))
main <- ""
else if (is.logical(main))
main <- if (main)
deparse(model$call)
else ""
else if (!is.character(main))
stop("`main' must be TRUE, FALSE, NULL or character
(vector).")
main <- rep(main, length = n.tms)
pf <- parent.frame()
carrier <- function(term) {
if (length(term) > 1)
carrier(term[[2]])
else eval(term, data, enclos = pf)
}
carrier.name <- function(term) {
if (length(term) > 1)
carrier.name(term[[2]])
else as.character(term)
}
if (is.null(xlabs))
xlabs <- unlist(lapply(cn, carrier.name))
if (partial.resid)
pres <- residuals(model, "partial")
is.fac <- sapply(nmt, function(i) is.factor(mf[, i]))
se.lines <- function(x, iy, i, ff = 2) {
tt <- ff * terms$se.fit[iy, i]
lines(x, tms[iy, i] + tt, lty = lty.se, lwd = lwd.se,
col = col.se)
lines(x, tms[iy, i] - tt, lty = lty.se, lwd = lwd.se,
col = col.se)
}
nb.fig <- prod(par("mfcol"))
if (ask) {
op <- par(ask = TRUE)
on.exit(par(op))
}
for (i in 1:n.tms) {
ylims <- range(tms[, i], na.rm = TRUE)
if (se)
ylims <- range(ylims, tms[, i] + 1.05 * 2 * terms$se.fit[,
i], tms[, i] - 1.05 * 2 * terms$se.fit[, i],
na.rm = TRUE)
if (partial.resid)
ylims <- range(ylims, pres[, i], na.rm = TRUE)
if (rug)
ylims[1] <- ylims[1] - 0.07 * diff(ylims)
if (is.fac[i]) {
ff <- mf[, nmt[i]]
ll <- levels(ff)
xlims <- range(seq(along = ll)) + c(-0.5, 0.5)
xx <- codes(ff)
if (rug) {
xlims[1] <- xlims[1] - 0.07 * diff(xlims)
xlims[2] <- xlims[2] + 0.03 * diff(xlims)
}
## use factor levels, not numbers, on x-axis
plot(1, 0, type = "n", xlab = xlabs[i], ylab = ylabs[i],
xlim = xlims, ylim = ylims, main = main[i],
xaxt="n",...)
axis(1,1:length(ll),ll)
for (j in seq(along = ll)) {
ww <- which(ff == ll[j])[c(1, 1)]
jf <- j + c(-0.4, 0.4)
lines(jf, tms[ww, i], col = col.term, lwd = lwd.term,
...)
if (se)
se.lines(jf, iy = ww, i = i)
}
}
else {
xx <- carrier(cn[[i]])
xlims <- range(xx, na.rm = TRUE)
if (rug)
xlims[1] <- xlims[1] - 0.07 * diff(xlims)
oo <- order(xx)
plot(xx[oo], tms[oo, i], type = "l", xlab = xlabs[i],
ylab = ylabs[i], xlim = xlims, ylim = ylims,
main = main[i], col = col.term, lwd = lwd.term,
...)
if (se)
se.lines(xx[oo], iy = oo, i = i)
}
if (partial.resid)
points(xx, pres[, i], cex = cex, pch = pch, col = col.res)
if (rug) {
n <- length(xx)
lines(rep(jitter(xx), rep(3, n)), rep(ylims[1] +
c(0, 0.05, NA) * diff(ylims), n))
if (partial.resid)
lines(rep(xlims[1] + c(0, 0.05, NA) * diff(xlims),
n), rep(pres[, i], rep(3, n)))
}
}
invisible(n.tms)
}
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at
stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._