[R] plot(cox.zph()): customize xlab & ylab
Adaikalavan Ramasamy
ramasamy at cancer.org.uk
Mon Jul 11 18:25:40 CEST 2005
Dan, I think this works fine now.
I think the 'ylab' argument in other plot function take only a single
value whereas in the case of plot.cox.zph, it needs a vector. It is
especially confusing when this function does plots everything on a
single page by default and all you see is the last plot.
Thus, can I suggest that you check that the length of the ylab is the
same as the number of terms in the coxph call. To do this you could put
something along the following anywhere after the line "yy <- x$y"
n <- length(coefficients(fit))
if( length(ylab) != n){
warning( paste("'ylab' is of length", length(ylab),
"and was expecting a length of", n, "elements.") )
}
Regards, Adai
On Mon, 2005-07-11 at 16:04 +0100, Dan Bebber wrote:
> Dear all,
>
> I've modified the plot.cox.zph function to allow
> customized xlab and ylab (see below). Someone might
> like to confirm that it works.
>
> Thanks for all the assistance.
>
> Dan
> ___________________________________
>
> plot.cox.zph <- function (x, resid = TRUE, se = TRUE,
> df = 4, nsmo = 40, var,
> xlab="Time",ylab = paste("Beta(t) for",
> dimnames(yy)[[2]]),...)
> {
> xx <- x$x
> yy <- x$y
> d <- nrow(yy)
> df <- max(df)
> nvar <- ncol(yy)
> pred.x <- seq(from = min(xx), to = max(xx), length
> = nsmo)
> temp <- c(pred.x, xx)
> lmat <- ns(temp, df = df, intercept = TRUE)
> pmat <- lmat[1:nsmo, ]
> xmat <- lmat[-(1:nsmo), ]
> qmat <- qr(xmat)
> if (se) {
> bk <- backsolve(qmat$qr[1:df, 1:df], diag(df))
> xtx <- bk %*% t(bk)
> seval <- d * ((pmat %*% xtx) * pmat) %*%
> rep(1, df)
> }
> if (missing(var))
> var <- 1:nvar
> else {
> if (is.character(var))
> var <- match(var, dimnames(yy)[[2]])
> if (any(is.na(var)) || max(var) > nvar ||
> min(var) <
> 1)
> stop("Invalid variable requested")
> }
> if (x$transform == "log") {
> xx <- exp(xx)
> pred.x <- exp(pred.x)
> }
> else if (x$transform != "identity") {
> xtime <- as.numeric(dimnames(yy)[[1]])
> apr1 <- approx(xx, xtime, seq(min(xx),
> max(xx), length = 17)[2 *
> (1:8)])
> temp <- signif(apr1$y, 2)
> apr2 <- approx(xtime, xx, temp)
> xaxisval <- apr2$y
> xaxislab <- rep("", 8)
> for (i in 1:8) xaxislab[i] <- format(temp[i])
> }
> for (i in var) {
> y <- yy[, i]
> yhat <- pmat %*% qr.coef(qmat, y)
> if (resid)
> yr <- range(yhat, y)
> else yr <- range(yhat)
> if (se) {
> temp <- 2 * sqrt(x$var[i, i] * seval)
> yup <- yhat + temp
> ylow <- yhat - temp
> yr <- range(yr, yup, ylow)
> }
> if (x$transform == "identity")
> plot(range(xx), yr, type = "n", xlab =
> xlab, ylab = ylab[i],
> ...)
> else if (x$transform == "log")
> plot(range(xx), yr, type = "n", xlab =
> xlab, ylab = ylab[i],
> log = "x", ...)
> else {
> plot(range(xx), yr, type = "n", xlab =
> xlab, ylab = ylab[i],
> axes = FALSE, ...)
> axis(1, xaxisval, xaxislab)
> axis(2)
> box()
> }
> if (resid)
> points(xx, y)
> lines(pred.x, yhat)
> if (se) {
> lines(pred.x, yup, lty = 2)
> lines(pred.x, ylow, lty = 2)
> }
> }
> }
>
>
> --- Thomas Lumley <tlumley at u.washington.edu> wrote:
>
> > On Mon, 11 Jul 2005, Adaikalavan Ramasamy wrote:
> >
> > > I am not sure if there is an easy way around this.
> > An ugly hack is to
> > > make a copy the function "survival:::plot.cox.zph"
> > and make your
> > > modified function. But there are others in the
> > list who might know
> > > neater solutions.
> >
> > If you then send a patch to the package maintainer
> > it stops being an ugly
> > hack and turns into an example of collaborative
> > open-source development :)
> >
> > -thomas
> >
> >
>
>
>
>
>
>
> ___________________________________________________________
> Yahoo! Messenger - NEW crystal clear PC to PC calling worldwide with voicemail http://uk.messenger.yahoo.com
>
More information about the R-help
mailing list