On Feb 3, 2011, at 3:27 PM, Bond, Stephen wrote:
Responding to the suggestion by D Winsemius :
s1 <-
survfit(mod1,newdata=inc[50050:50100,c("strt","stp","incpost",
"amt","rate", "termfac")],
+ se.fit=F,individual=T,type="aa")
Error in Surv(time = strt, time2 = stp, event = (resp == 1)) :
object 'resp' not found
it appears it wants to fit a survival curve instead of predicting a
survival curve for the subject using the explanatory vars. I think I
probably need
predict.coxph(...,type="expected")
or write code for matching time dependent risk to the proper time
index of the survival curve.
I would have interpreted that error message as saying the function
also needed some of the information on the LHS of the formula. Keep
feeding it data that it complains is missing.
Thanks David.
Stephen
-----Original Message-----
From: David Winsemius [mailto:dwinsem...@comcast.net]
Sent: Thursday, February 03, 2011 3:10 PM
To: Bond, Stephen
Cc: r-help@r-project.org
Subject: Re: [R] coxph fails to survfit
On Feb 3, 2011, at 2:14 PM, Bond, Stephen wrote:
I have a model with quant vars only and the error message does not
make sense:
(mod1 <- coxph(Surv(time=strt,time2=stp,event=(resp==1))~ +incpost
+I(amt/1e5)+rate+strata(termfac),
subset=dt<"2010-08-30", data=inc,method="efron"))
Call:
coxph(formula = Surv(time = strt, time2 = stp, event = (resp ==
1)) ~ +incpost + I(amt/1e+05) + rate + strata(termfac), data = inc,
subset = dt < "2010-08-30", method = "efron")
coef exp(coef) se(coef) z p
incpost 0.2563 1.292 0.02479 10.34 0.0e+00
I(amt/1e+05) -0.0532 0.948 0.00487 -10.92 0.0e+00
rate -0.0507 0.951 0.00945 -5.36 8.2e-08
Likelihood ratio test=295 on 3 df, p=0 n= 1192634
length(mod1$xlevels)
[1] 0
# now calling survfit with just a few rows from the data
s1 <-
survfit
(mod1,newdata=inc[50050:50100,],se.fit=F,individual=T,type="aa")
Error in `contrasts<-`(`*tmp*`, value = "contr.treatment") :
contrasts can be applied only to factors with 2 or more levels
I obviously cannot test my theory, but would have thought the call
would be:
s1 <- survfit(mod1, newdata=inc[50050:50100,], c("incpost", "amt,
"rate", "termfac") ],
se.fit=F,individual=T,type="aa")
I.e., a data.frame for newdata that matched the variables used in the
RHS of the formula.
inc[50050:50100,]
dt inc incpost strt stp resp rate amt p1 matfac
termfac vol
50050 2009-02-05 0.00000 0.00 0 1 0 5.35 266833 C
5 3 14229571
50051 2009-02-06 -0.09575 0.00 1 2 0 5.35 266833 C
5 3 14229571
50052 2009-02-09 -0.03875 0.00 4 5 0 5.35 266833 C
5 3 14229571
50053 2009-02-10 0.00850 0.00 5 6 0 5.35 266833 C
5 3 14229571
50054 2009-02-11 0.00250 0.00 6 7 0 5.35 266833 C
5 3 14229571
50055 2009-02-12 -0.04725 0.00 7 8 0 5.35 266833 C
5 3 14229571
50056 2009-02-13 -0.00075 0.00 8 9 0 5.35 266833 C
5 3 14229571
50057 2009-02-16 -0.10425 0.00 11 12 0 5.35 266833 C
5 3 14229571
50058 2009-02-17 -0.10425 0.00 12 13 0 5.35 266833 C
5 3 14229571
50059 2009-02-18 -0.02925 0.00 13 14 0 5.35 266833 C
5 3 14229571
50060 2009-02-19 -0.01300 0.00 14 15 0 5.35 266833 C
5 3 14229571
50061 2009-02-20 -0.04075 0.00 15 16 0 5.35 266833 C
5 3 14229571
50062 2009-02-23 -0.03100 0.00 18 19 0 5.35 266833 C
5 3 14229571
50063 2009-02-24 0.01975 0.00 19 20 0 5.35 266833 C
5 3 14229571
50064 2009-02-25 0.13050 0.00 20 21 0 5.35 266833 C
5 3 14229571
50065 2009-02-26 0.12725 0.00 21 22 0 5.35 266833 C
5 3 14229571
... further output truncated
Please, comment if you see what's the issue.
Thank you.
Stephen
[[alternative HTML version deleted]]
______________________________________________
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.
David Winsemius, MD
West Hartford, CT
David Winsemius, MD
West Hartford, CT
______________________________________________
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.
