Dear all,
I was wondering if anyone can help? I am an R user but recently I have resorted
to SAS to calculate the probability of the event (and the associated confidence
interval) for the Cox model with combinations of risk factors. For example,
suppose I have a Cox model with two binary variables, one for gender and one for
treatment, I wish to calculate the probability of survival for the combinations
of 1s and 0s for both variables. In SAS the following code would be used.
I have used the code below to obtain the survival estimates (at 1 year) and
confidence intervals for combinations of risk factors for my outcome.
/* Combinations of Risk Factors */
data test2;
input sex treat;
DATALINES;
0 0
1 0
0 1
1 1
;
run;
/* Survival estimates for the above combinations */
proc phreg data = pudat2;
model withtime*wcens(0) = sex treat /ties = efron;
baseline out = surv2 survival = survival lower = slower upper = supper
covariates = test2 /method = ch nomean cltype=loglog;
run;
/* Survival estimates at 1 year */
proc print data = surv2 noobs;
where withtime = 364;
run;
I now would like to do the same thing but in a competing risks setting
(cumulative incidence rather than Cox model version (crr from cmprsk in R)
however there is no pre-set routine for this in SAS. I am therefore keen to use
R to obtain the survival estimates for both the overall outcome and the
competing risks outcomes.
Can anyone help me to translate the SAS code above into R code so that I can
check I am doing things correctly? I am using R version 2.9.2 on Windows XP
Professional.
I think it would begin as follows:
library(cmprsk)
fit1 <- coxph(Surv(withtime,wcens)~ sex+treat,data=pudat2)
my.fit <- survift(fit1)
summary(my.fit)
Unfortunately, the estimates shown do not match those of SAS.
Additionally, can you use survfit with crr? If not, any suggestions as to how I
may obtain the necessary estimates in a competing risks setting?
Thank you for any help you can give.
Laura
[[alternative HTML version deleted]]
-- Begin included message
/* Combinations of Risk Factors */
data test2;
input sex treat;
DATALINES;
0 0
1 0
0 1
1 1
;
run;
/* Survival estimates for the above combinations */
proc phreg data = pudat2;
model withtime*wcens(0) = sex treat /ties = efron;
baseline out = surv2 survival = survival lower = slower upper = supper
covariates = test2 /method = ch nomean cltype=loglog;
run;
/* Survival estimates at 1 year */
proc print data = surv2 noobs;
where withtime = 364;
fit1 <- coxph(Surv(withtime,wcens)~ sex+treat,data=pudat2)
my.fit <- survift(fit1)
summary(my.fit)
Unfortunately, the estimates shown do not match those of SAS.
-- End inclusion
You need to add a newdata=test2 argument to the survfit call. Like the
SAS code, test2 is a data set containing the list of subjects for which
you want an estimate. You can also add a times= argument to the summary
function.
Last, if you look at the documentation you will see that the default
handling for ties in the baseline hazard is the Fleming-Harrington
method in R, whenever the "Efron" approx was used in computing the
coefficients. (Same approximation, 2 different names when used in the
two contexts, for more details see the chapter in my book). The "CH"
or
Aalen approx for a baseline hazard is the same math as the Breslow
approx for the coefficients. So to match the chimeric mixture you used
in the SAS code you will need to add type="aalen" to the survfit
call.
The numeric effect of this is usually neglible, but when I see "do not
match" in a query I assume you are looking to have all the digits the
same.
As to competing risks, I have some code for R but documentation is yet
to be done -- on my to do list. Walter Kremers
(kremers.walter at mayo.edu) has a nice SAS macro for this case, however.
You might want to drop him a line.
Terry Therneau