Dimitri Szerman
2006-Jun-12 22:00 UTC
[R] non parametric estimates of the hazard with right censored data
Hi, I want to plot non parametric estimates of the empirical hazard function for right censored data. I've tried many functions from different packages (muhaz, Design, survival, eha, event), but none of them gave me what I wanted. Am I missing something? Here's what I want. The data below is the same used by Kiefer (J. Economic Literature, 1988), which in turns use a subset of the data used by Kiefer (J. Econometrics, 1985) in his study of strikes duration. ## R dur = c( 7, 9, 13 ,14 ,26 ,29 ,52 ,80 ,9 ,37 ,41 ,49 ,52 ,80 ,3 ,17 ,19 ,28 ,72 ,80 ,80 ,80 ,80 ,80 ,80 ,15 ,61,80,2,25 ,80 ,3 ,10 ,1 ,2 ,2 ,3 ,3 ,4 ,8,11 ,22 ,23 ,27 ,32 ,33 ,35 ,43 ,43 ,44 ,80 ,5 ,49 ,2 ,12 ,12 ,21 ,21 ,27 ,38 ,42 ,80) cen = c(1, 1 ,1 ,1 ,1 ,1 ,1 ,0,1 ,1 ,1 ,1 ,1 ,0 ,1 ,1 ,1 ,1 ,1 ,0 ,0 ,0 ,0 ,0 ,0 ,1 ,1 ,0 ,1 ,1 ,0 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,0 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1,1, 0) library(survival) hfit = survfit( Surv(dur, cen) ~ 1 ) tevent = data.frame( tj=summary(hfit)$time, hj=summary(hfit)$n.event, nj=summary(hfit)$n.risk, haz=summary(hfit)$n.event/summary(hfit)$n.risk) plot(tevent$tj, tevent$haz, type="s", xlim=c(0,80), ylim=c(0,.13) ) ## END This plot is the same presented by Kiefer (figure 7). I even managed to plot the *cumulative*, or *integrated* hazard (figure 8), but the hazard itself. As you can see, I could do it myself, but as I want to do (many, many) plots by strata, a ready-to-go function would help a lot. Thany you so much, Dimitri [[alternative HTML version deleted]]