Xinyi Li asked how to keep track of which coxph models, called from within a
loop, were responsible for warning messges. One solution is to modify the
coxph code so that those models are marked in the return coxph object. Below
is a set of changes to the final 40 lines of coxph.fit, that will cause the
component "infinite.warn" to be added to the result whenever a warning was
generated; it will be a vector of T/F showing which component(s) of the
coefficient vector generated the warning.
Change the code for coxph.fit.s, then do
> source('coxph.fit.revised.s') #or whatever you called it
> coxph <- source('coxph.s')
> coxph.wtest <- survival:::coxph.wtest
Line 2 causes you to have a local version of coxph, otherwise, due
to name spaces, the original version of coxph.fit will still be used. Line 3
is another consequence of name spaces. Then code such as
for(i in 1:ncol(x)){
fit=coxph(TIME~x[,i])
if (!is.null(fit$infinite.warn)) cat("Waring in fit", i, "\n")
sfit <- summary(fit)
results[i,1]=sfit$coef[1]
results[i,2]=sfit$coef[3]
results[i,3]=sfit$coef[5]
}
Will report out the models with a warning.
Notes
1. The warning message is not completely reliable
2. Name spaces protect a package from accidental overrides, when a user or
some other package reuses a function name. With its hundreds of packages, they
are a necessity for R. But sometimes you really do want to override a
function; then they are a bit of a pain. Others with a better grasp of R
internals may be able to suggest a better override than I have done here.
Terry Therneau
-------------------------
infs <- abs(coxfit$u %*% var)
keep.infs <- F # new line
if (maxiter >1) {
if (coxfit$flag == 1000)
warning("Ran out of iterations and did not converge")
else {
infs <- ((infs > control$eps) &
infs > control$toler.inf*abs(coef))
if (any(infs))
warning(paste("Loglik converged before variable ",
paste((1:nvar)[infs],collapse=","),
"; beta may be infinite. "))
keep.infs <- T #new line
}
}
names(coef) <- dimnames(x)[[2]]
lp <- c(x %*% coef) + offset - sum(coef*coxfit$means)
score <- exp(lp[sorted])
coxres <- .C("coxmart", as.integer(n),
as.integer(method=='efron'),
stime,
sstat,
newstrat,
as.double(score),
as.double(weights),
resid=double(n))
resid <- double(n)
resid[sorted] <- coxres$resid
names(resid) <- rownames
coef[which.sing] <- NA
temp <- list(coefficients = coef, #modified line
var = var,
loglik = coxfit$loglik,
score = coxfit$sctest,
iter = coxfit$iter,
linear.predictors = as.vector(lp),
residuals = resid,
means = coxfit$means,
method='coxph')
if (keep.infs) temp$infinite.warn <- infs #new line
temp #new line
}
}
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
