[R] Problem with survfit when model has interaction with strata
Sebastian Pölsterl
sebp at k-d-w.org
Thu Aug 6 15:44:38 CEST 2015
Dear all,
I encountered a problem when running survfit using a coxph model that
contains an interaction with a strata variable (see the attached
example). I hope this is the right place to report problems.
> source("~/survfit-example.R")
Error in scale.default(x2, center = xcenter, scale = FALSE) :
length of 'center' must equal the number of columns of 'x'
> traceback()
9: stop("length of 'center' must equal the number of columns of 'x'")
8: scale.default(x2, center = xcenter, scale = FALSE)
7: scale(x2, center = xcenter, scale = FALSE)
6: survfit.coxph(fit, newdata = newd)
5: survfit(fit, newdata = newd) at survfit-example.R#10
4: eval(expr, envir, enclos)
3: eval(ei, envir)
2: withVisible(eval(ei, envir))
1: source("~/survfit-example.R")
> sessionInfo()
R version 3.2.1 (2015-06-18)
Platform: x86_64-unknown-linux-gnu (64-bit)
Running under: Fedora 22 (Twenty Two)
locale:
[1] C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] survival_2.38-3 foreign_0.8-65
loaded via a namespace (and not attached):
[1] tools_3.2.1 splines_3.2.1
I believe the problem is occurring on line 263 of survfit.coxph.R, where
a model.matrix from newdata is created. The values passed to scale look
like:
> x2
RXplacebo RXtreatment:strata(SEX)male RXplacebo:strata(SEX)male
1 0 1 0
2 0 0 0
> xcenter
RXplacebo strata(SEX)male:RXplacebo
0.5000000 0.2380952
Clearly, model.matrix is adding an additional (redundant) column that
was omitted when fitting the original model. Since interactions with
strata variables are a special cases in coxph, I believe similar care
should be taken when creating the model matrix for newdata.
Best regards,
Sebastian Pölsterl
-------------- next part --------------
library(foreign)
library(survival)
data <- data.frame(read.spss("http://web1.sph.emory.edu/dkleinb/allDatasets/surv2datasets/anderson.sav"))
S <- with(data, Surv(SURVT, I(STATUS == "relapse")))
fit <- coxph(S ~ strata(SEX)*RX, data=data)
newd <- data.frame(RX=c("treatment", "treatment"), SEX=c("male", "female"))
sfit <- survfit(fit, newdata=newd)
More information about the R-help
mailing list