Search code examples
rcox-regressionsurvival

Cox proportional hazard model on (burn) data


On a study is described which evaluates a protocol change in disinfectant practices in a large midwestern university medical center. Of primary interest in the study is a comparison of two methods of body cleansing. The first method, used exclusively from January 1983 to June 1984, consisted of a routine bathing care method (initial surface decontamination with 10% povidone-iodine followed with regular bathing with Dial soap). From June 1984 to the end of the study period in December1985, body cleansing was initially performedusing 4% chlorhexidine gluconate. Eighty-four patients were in the group who received the new bathing solution, chlorhexidine, and 70 patients served as the control group who received routine bathing care, povidoneiodine. Included in the data set is a covariate that measures the total surface area burned. The data is (burn). I want to test for: 1- any difference in survival functions for the two groups. 2- any difference in survival functions for the two groups adjusting for total area burned.

library(KMsurv)
data()
data(burn)
burn

library(survival)

I know that the function that would be used is coxph(), but I'm not sure which groups that I should test for (from the above information). Are they T1 and D2? so that for 1, Coxfit1<-coxph(Surv(T1,D2)~group, data = burn)? and for 2, Coxfit2<-coxph(Surv(T1,D2)~Z4, data = burn)?

What is this code doing?

for(i in 1:154){
  if (burn$??[i]==2)
    burn$Z1[i]<-1
  else burn$Z1[i]<-0
}

for(i in 1:154){
  if (burn$??[i]==3)
    burn$Z2[i]<-1
  else burn$Z2[i]<-0
}


Solution

  • For question 1, you want to test the survival distributions* between the levels of the Z1 variable. There is no variable called group in the dataset. Z1=0 means routine bathing and Z1=1 means body cleansing. You may want to convert all Z variables to factors before continuing further (except Z4).

    library(survival)
    library(KMsurv)
    library (dplyr)
    
    burn$Z1 <- factor(burn$Z1, label=c("Routine bathing", "Body cleansing"))
    

    * The word survival needs some clarification. Presumably it is time until first straphylocous aureaus infection (D3) or on study time if no event occurred. The time is in variable T3.

    The command to perform the test is:

    coxph(Surv(T3,D3) ~ Z1, data=burn)
                        coef exp(coef) se(coef)      z      p
    Z1Body cleansing -0.5614    0.5704   0.2934 -1.914 0.0557
    

    For question 2, Z4 contains the percentage of total surface area burned, the variable to adjust for.

    coxph(Surv(T3,D3)~Z1+Z4, data=burn)
    
                          coef exp(coef)  se(coef)      z     p
    Z1Body cleansing -0.524764  0.591695  0.295769 -1.774 0.076
    Z4                0.007248  1.007275  0.007145  1.015 0.310
    

    So there appears to be no difference in time until first infection between those who were given routine bathing vs body cleansing.