[R] coxph: cumulative mortality hazard over time with associated confidence intervals

[Lise Aubry]

Dear R-users,

I computed a simple coxph model and plotted survival over time with
associated confidence intervals for 2 covariate levels (males and
females).

M1  <- coxph(survobject~sex, data=surv)
M1
survsex <- survfit(survobject~sex,data=surv)
summary(survsex)
plot(survsex, conf.int=T, col=c("black","red"), lty = c(1,2),
lwd=c(1,2), xlab="Time", ylab="Survival probability")
legend(10, 0.2, c("males","females"), col = c("black","red"), lty =
c(1,2), lwd=c(1,2),bty='n')

I would like to make a similar plot of the cumulative mortality hazard
over time for both males and females with associated confidence
intervals. I would like to make sure I obtained proper estimates of
the cumulative hazard for both males and females with associated
confidence intervals but I have doubts regarding my calculations

# I used basehaz ‘H0’ to get to the baseline hazard (covariate level
0, males in this case)
basehazard <- basehaz(M1, centered=FALSE)[,1]
 NBLOWCI <- basehazard * exp(1.96*sqrt(M1$var[1]))
 NBHIGHCI <- basehazard *  exp(-1.96*sqrt(M1$var[1]))

# I applied the following equation H0*exp(Beta*X) to get to the
cumulative hazard over time (covariate level 1) for females, In # this
case X=1, and Beta is the regression coefficient from the coxph model
(M1$coef[1])
 bandhazard <- basehazard * exp(M1$coef[1]*1)
 BLOWCI <- bandhazard * exp(1.96*sqrt(M1$var[1]))
 BHIGHCI <- bandhazard *  exp(-1.96*sqrt(M1$var[1]))

# then the cumulative hazard plot with associated confidence intervals
for males and females
 Time <-  basehaz(M1, centered=FALSE)[,2]
 plot(Time, basehazard, ylim=c(0,3), type='l', col="black", lty=1,
lwd=1, xlab="Time”, ylab="Mortality  Hazard")
 lines(Time,NBLOWCI,type='l', col="black", lty=2, lwd=1)
  lines(Time,NBHIGHCI,type='l', col="black", lty=2, lwd=1)
  lines(Time, bandhazard, col="red", lty=1, lwd=1)
  lines(Time,BLOWCI,type='l', col="red", lty=2, lwd=1)
  lines(Time,BHIGHCI,type='l', col="red", lty=2, lwd=1)
  legend(15, 4, c("Males","Females"), col = c("black","red"), lty =
c(1,2), lwd=c(1,2),bty='n')

1.	Is this correct?
2.	Is there any automated way to do this in R?
3.	 Similarly, is there any automated way to plot the mortality hazard
(instead of the cumulative hazard) over time with associated
confidence intervals for males and females?

Any comments or suggestions are welcome!

Lise