[R] coxph() - unexpected result using Crawley's seedlings data (The R Book)
Jacob Brogren
jacob at brogren.nu
Tue Jun 28 17:04:20 CEST 2011
Hi,
sorry about that; here is the full output - data set, structure, model and result.
Cheers
Jacob
> seedlings
cohort death gapsize status
1 September 7 0.5889 1
2 September 3 0.6869 1
3 September 12 0.1397 1
4 September 1 0.1921 1
5 September 4 0.2798 1
6 September 2 0.2607 1
7 September 6 0.9467 1
8 September 6 0.6375 1
9 September 8 0.1527 1
10 September 3 0.8237 1
11 September 1 0.5979 1
12 September 1 0.2914 1
13 September 3 0.5053 1
14 September 5 0.4714 1
15 September 2 0.6041 1
16 September 8 0.8812 1
17 September 4 0.8416 1
18 September 1 0.0577 1
19 September 2 0.9034 1
20 September 2 0.4348 1
21 September 7 0.9878 1
22 September 11 0.1486 1
23 September 14 0.5003 1
24 September 1 0.8507 1
25 September 10 0.8187 1
26 September 14 0.0291 1
27 September 1 0.3785 1
28 September 4 0.8384 1
29 September 2 0.8351 1
30 September 2 0.9674 1
31 October 1 0.6943 1
32 October 1 0.2591 1
33 October 2 0.7397 1
34 October 2 0.4663 1
35 October 14 0.9115 1
36 October 5 0.1750 1
37 October 1 0.5628 1
38 October 8 0.2681 1
39 October 5 0.6967 1
40 October 2 0.7020 1
41 October 4 0.7971 1
42 October 3 0.4047 1
43 October 5 0.0498 1
44 October 10 0.0364 1
45 October 9 0.4080 1
46 October 1 0.6226 1
47 October 11 0.3002 1
48 October 3 0.8111 1
49 October 21 0.4894 1
50 October 1 0.0375 1
51 October 4 0.2560 1
52 October 9 0.2168 1
53 October 8 0.7437 1
54 October 1 0.9082 1
55 October 3 0.9496 1
56 October 9 0.1040 1
57 October 9 0.8691 1
58 October 16 0.9502 1
59 October 6 0.0790 1
60 October 1 0.5658 1
> str(seedlings)
'data.frame': 60 obs. of 4 variables:
$ cohort : Factor w/ 2 levels "October","September": 2 2 2 2 2 2 2 2 2 2 ...
$ death : int 7 3 12 1 4 2 6 6 8 3 ...
$ gapsize: num 0.589 0.687 0.14 0.192 0.28 ...
$ status : num 1 1 1 1 1 1 1 1 1 1 ...
> model1 <- coxph(Surv(death,status)~strata(cohort)*gapsize,data=seedlings)
> 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
28 jun 2011 kl. 15.48 skrev Robert A LaBudde:
> 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