I'd like to calculate expected risk (cumulative incidences), which are derived from fitted Cox PH model using R packages.
I have the fitted Cox PH model like as follows:
[Variables] Dataset: 10,000 cases(patients) with 6 variabels
t: time to event, s: event (coded as 0 or 1)
covariates: X1, X2, X3, X4 (coded as 0 or 1)
fit <- coxph (Surv (t,s) ~ X1+X2+X3+X4, data=data)
coef exp(coef) se (coef) z P X1 -0.3777 0.6855 0.1120 -3.37 0.00075 X2 0.4014 1.4938 0.0518 7.74 <0.0001 X3 0.7417 2.0995 0.0893 8.31 <0.0001 X4 0.4330 1.5419 0.1268 3.42 <0.001
From this model, I'd like to calculate expected risk (cumulative incidence of events) for each cases according to X1 = 0 or X1 = 1.
In other words, if the X1 of all 10,000 cases were 0, how could I calculate expected risk for each cases? At the same time, if the X1 of all 10,000 cases were 1, how could I calculate expected risk for each cases? (using R)
After calculating the expected risk for each patients, then I'd like to calculate the risk-benefit ratio for each patients according to variable X1 [by (Expected risk, when X1 = 0) - (Expected risk, when X1 = 1)], dose it right?
Meanwhile, I have tried to plot the expected cumulative risk curve according to X1 like as follows; was it appropriate ?
baseha = basehaz (fit, centered=FALSE)
X1=0: exp (-0.38*0) = 1, X1=1: exp (-0.38*1) = 0.68
plot (baseha$$time[,2], baseha$hazard*(1), type="s", lty=1)
lines (baseha$$time[,2], baseha$hazard*(0.68), type="s", lty=2)
Thanks for your kind answers.