[R] coxph() - unexpected result using Crawley's seedlings data (The R Book)

Robert A LaBudde ral at lcfltd.com
Tue Jun 28 15:48:32 CEST 2011


Did you create the 'status' variable the way indicated on p. 797?

Frequently with Surv() it pays to use syntax such as Surv(death, 
status==1) to make a clear logical statement of what is an event 
(status==1) vs. censored.

PS. Next time include head(seedlings) and str(seedlings) to make 
clear what you are using as data.


At 06:51 AM 6/28/2011, Jacob Brogren wrote:
>Hi,
>
>I ran the example on pp. 799-800 from Machael Crawley's "The R Book" 
>using package survival v. 2.36-5, R 2.13.0 and RStudio 0.94.83. The 
>model is a Cox's Proportional Hazards model. The result was quite 
>different compared to the R Book. I have compared my code to the 
>code in the book but can not find any differences in the function 
>call. My results are attached as well as a link to the results 
>presented in the book (link to Google Books).
>
>When running the examples on pp. 797-799 I can't detect any 
>differences in results so I don't think there are errors in the data 
>set or in the creation of the status variable.
>
>---------------------------------
>Original from the R Book:
>http://books.google.com/books?id=8D4HVx0apZQC&lpg=PA799&ots=rQgd_8ofeS&dq=r%20coxph%20crawley&pg=PA799#v=onepage&q&f=false
>
>---------------------------------
>My result:
> > summary(model1)
>Call:
>coxph(formula = Surv(death, status) ~ strata(cohort) * gapsize,
>     data = seedlings)
>
>   n= 60, number of events= 60
>
>                                             coef 
> exp(coef)  se(coef)      z Pr(>|z|)
>gapsize                                -0.001893  0.998109  0.593372 
>-0.003    0.997
>gapsize:strata(cohort)cohort=September  0.717407  2.049112  0.860807 
>  0.833    0.405
>
>                                        exp(coef) exp(-coef) lower 
> .95 upper .95
>gapsize                                   0.9981      1.002 
>0.3120     3.193
>gapsize:strata(cohort)cohort=September    2.0491      0.488 
>0.3792    11.074
>
>Rsquare= 0.022   (max possible= 0.993 )
>Likelihood ratio test= 1.35  on 2 df,   p=0.5097
>Wald test            = 1.32  on 2 df,   p=0.5178
>Score (logrank) test = 1.33  on 2 df,   p=0.514
>
>Anyone have an idea why this is occurring?
>
>Kind Regards
>
>Jacob
>______________________________________________
>R-help at r-project.org mailing list
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>and provide commented, minimal, self-contained, reproducible code.

================================================================
Robert A. LaBudde, PhD, PAS, Dpl. ACAFS  e-mail: ral at lcfltd.com
Least Cost Formulations, Ltd.            URL: http://lcfltd.com/
824 Timberlake Drive                     Tel: 757-467-0954
Virginia Beach, VA 23464-3239            Fax: 757-467-2947

"Vere scire est per causas scire"



More information about the R-help mailing list