R Under development (unstable) (2023-05-10 r84417) -- "Unsuffered Consequences"
Copyright (C) 2023 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> options(na.action=na.exclude) # preserve missings
> options(contrasts=c('contr.treatment', 'contr.poly')) #ensure constrast type
> library(survival)
>
> expect <- survexp(futime ~ 1,
+ rmap = list(age=(accept.dt - birth.dt), sex=1,
+ year=accept.dt, race='white'), jasa, cohort=F,
+ ratetable=survexp.usr)
>
> survdiff(Surv(jasa$futime, jasa$fustat) ~ offset(expect))
Call:
survdiff(formula = Surv(jasa$futime, jasa$fustat) ~ offset(expect))
Observed Expected Z p
75.000 0.587 -97.119 0.000
> # Now fit the 6 models found in Kalbfleisch and Prentice, p139
> sfit.1 <- coxph(Surv(start, stop, event)~ (age + surgery)*transplant,
+ jasa1, method='breslow')
> sfit.2 <- coxph(Surv(start, stop, event)~ year*transplant,
+ jasa1, method='breslow')
> sfit.3 <- coxph(Surv(start, stop, event)~ (age + year)*transplant,
+ jasa1, method='breslow')
> sfit.4 <- coxph(Surv(start, stop, event)~ (year +surgery) *transplant,
+ jasa1, method='breslow')
> sfit.5 <- coxph(Surv(start, stop, event)~ (age + surgery)*transplant + year ,
+ jasa1, method='breslow')
> sfit.6 <- coxph(Surv(start, stop, event)~ age*transplant + surgery + year,
+ jasa1, method='breslow')
>
> summary(sfit.1)
Call:
coxph(formula = Surv(start, stop, event) ~ (age + surgery) *
transplant, data = jasa1, method = "breslow")
n= 170, number of events= 75
coef exp(coef) se(coef) z Pr(>|z|)
age 0.01386 1.01395 0.01813 0.765 0.445
surgery -0.54652 0.57896 0.61091 -0.895 0.371
transplant 0.11572 1.12268 0.32729 0.354 0.724
age:transplant 0.03473 1.03534 0.02725 1.274 0.202
surgery:transplant -0.29037 0.74799 0.75819 -0.383 0.702
exp(coef) exp(-coef) lower .95 upper .95
age 1.014 0.9862 0.9786 1.051
surgery 0.579 1.7272 0.1748 1.917
transplant 1.123 0.8907 0.5911 2.132
age:transplant 1.035 0.9659 0.9815 1.092
surgery:transplant 0.748 1.3369 0.1692 3.306
Concordance= 0.595 (se = 0.037 )
Likelihood ratio test= 12.45 on 5 df, p=0.03
Wald test = 11.62 on 5 df, p=0.04
Score (logrank) test = 12.02 on 5 df, p=0.03
> sfit.2
Call:
coxph(formula = Surv(start, stop, event) ~ year * transplant,
data = jasa1, method = "breslow")
coef exp(coef) se(coef) z p
year -0.2652 0.7671 0.1051 -2.522 0.0117
transplant -0.2871 0.7504 0.5137 -0.559 0.5762
year:transplant 0.1371 1.1469 0.1409 0.973 0.3306
Likelihood ratio test=8.61 on 3 df, p=0.03488
n= 170, number of events= 75
> summary(sfit.3)
Call:
coxph(formula = Surv(start, stop, event) ~ (age + year) * transplant,
data = jasa1, method = "breslow")
n= 170, number of events= 75
coef exp(coef) se(coef) z Pr(>|z|)
age 0.01558 1.01571 0.01734 0.899 0.36887
year -0.27413 0.76023 0.10588 -2.589 0.00962 **
transplant -0.59388 0.55218 0.54222 -1.095 0.27339
age:transplant 0.03380 1.03438 0.02795 1.209 0.22653
year:transplant 0.20228 1.22419 0.14247 1.420 0.15566
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
exp(coef) exp(-coef) lower .95 upper .95
age 1.0157 0.9845 0.9818 1.0508
year 0.7602 1.3154 0.6178 0.9356
transplant 0.5522 1.8110 0.1908 1.5981
age:transplant 1.0344 0.9668 0.9792 1.0926
year:transplant 1.2242 0.8169 0.9259 1.6185
Concordance= 0.63 (se = 0.036 )
Likelihood ratio test= 14.85 on 5 df, p=0.01
Wald test = 13.77 on 5 df, p=0.02
Score (logrank) test = 14.06 on 5 df, p=0.02
> sfit.4
Call:
coxph(formula = Surv(start, stop, event) ~ (year + surgery) *
transplant, data = jasa1, method = "breslow")
coef exp(coef) se(coef) z p
year -0.2541 0.7756 0.1077 -2.360 0.0183
surgery -0.2365 0.7894 0.6282 -0.377 0.7065
transplant -0.2969 0.7431 0.5054 -0.588 0.5568
year:transplant 0.1653 1.1797 0.1416 1.167 0.2431
surgery:transplant -0.5499 0.5770 0.7759 -0.709 0.4785
Likelihood ratio test=12.36 on 5 df, p=0.03021
n= 170, number of events= 75
> sfit.5
Call:
coxph(formula = Surv(start, stop, event) ~ (age + surgery) *
transplant + year, data = jasa1, method = "breslow")
coef exp(coef) se(coef) z p
age 0.01504 1.01516 0.01760 0.855 0.3928
surgery -0.42025 0.65689 0.61564 -0.683 0.4948
transplant 0.07412 1.07694 0.33113 0.224 0.8229
year -0.13631 0.87258 0.07097 -1.921 0.0548
age:transplant 0.02691 1.02728 0.02712 0.992 0.3210
surgery:transplant -0.29656 0.74337 0.75800 -0.391 0.6956
Likelihood ratio test=16.2 on 6 df, p=0.0127
n= 170, number of events= 75
> sfit.6
Call:
coxph(formula = Surv(start, stop, event) ~ age * transplant +
surgery + year, data = jasa1, method = "breslow")
coef exp(coef) se(coef) z p
age 0.01527 1.01539 0.01750 0.873 0.3828
transplant 0.04457 1.04558 0.32171 0.139 0.8898
surgery -0.62113 0.53734 0.36786 -1.688 0.0913
year -0.13608 0.87278 0.07090 -1.919 0.0550
age:transplant 0.02703 1.02740 0.02714 0.996 0.3193
Likelihood ratio test=16.06 on 5 df, p=0.006687
n= 170, number of events= 75
>
> # Survival curve for an "average" subject,
> # done once as overall, once via individual method
> surv1 <- survfit(sfit.1, newdata=list(age=-2, surgery=0, transplant=0))
> newdata <- data.frame(start=c(0,50,100), stop=c(50,100, max(jasa1$stop)),
+ event=c(1,1,1), age=rep(-2,3), surgery=rep(0,3),
+ transplant=rep(0,3), name=c("Smith", "Smith", "Smith"))
> surv2 <- survfit(sfit.1, newdata, id=name)
> # Have to use unclass to avoid [.survfit trying to pick curves,
> # remove the final element "call" because it won't match, nor will newdata
> ii <- match(c("newdata", "call"), names(surv1))
> all.equal(unclass(surv1)[-ii],
+ unclass(surv2)[-ii])
[1] TRUE
>
>
> # Survival curve for a subject of age 50, with prior surgery, tx at 6 months
> # Remember that 'age' in jasa 1 was centered at 48
> data <- data.frame(start=c(0,183), stop=c(183,3*365), event=c(1,1),
+ age=c(2,2), surgery=c(1,1), transplant=c(0,1), id=c(1,1))
> summary(survfit(sfit.1, data, id=id))
Call: survfit(formula = sfit.1, newdata = data, id = id)
time n.risk n.event survival std.err lower 95% CI upper 95% CI
0.5 103 1 0.994 0.00722 0.980 1.000
1.0 102 3 0.975 0.01860 0.939 1.000
2.0 99 3 0.956 0.02914 0.900 1.000
4.0 96 2 0.943 0.03605 0.875 1.000
5.0 94 2 0.930 0.04286 0.849 1.000
7.0 92 1 0.923 0.04623 0.837 1.000
8.0 91 1 0.917 0.04959 0.824 1.000
11.0 89 1 0.910 0.05294 0.812 1.000
15.0 88 3 0.890 0.06278 0.775 1.000
16.0 85 1 0.883 0.06608 0.763 1.000
17.0 84 1 0.877 0.06928 0.751 1.000
20.0 83 2 0.864 0.07538 0.728 1.000
27.0 81 1 0.857 0.07849 0.716 1.000
29.0 80 1 0.850 0.08160 0.705 1.000
31.0 78 1 0.844 0.08473 0.693 1.000
34.0 77 1 0.837 0.08786 0.681 1.000
35.0 76 1 0.830 0.09098 0.669 1.000
36.0 75 1 0.823 0.09412 0.658 1.000
38.0 74 1 0.816 0.09727 0.646 1.000
39.0 72 2 0.802 0.10349 0.623 1.000
42.0 70 1 0.795 0.10664 0.611 1.000
44.0 69 1 0.788 0.10982 0.600 1.000
49.0 68 1 0.781 0.11300 0.588 1.000
50.0 67 1 0.774 0.11614 0.577 1.000
52.0 66 1 0.767 0.11925 0.565 1.000
57.0 65 1 0.760 0.12238 0.554 1.000
60.0 64 1 0.752 0.12552 0.542 1.000
65.0 63 1 0.745 0.12866 0.531 1.000
67.0 62 2 0.730 0.13494 0.508 1.000
68.0 60 1 0.722 0.13809 0.497 1.000
71.0 59 2 0.707 0.14420 0.474 1.000
76.0 57 1 0.699 0.14729 0.463 1.000
77.0 56 1 0.691 0.15043 0.451 1.000
79.0 55 1 0.683 0.15362 0.439 1.000
80.0 54 1 0.674 0.15680 0.428 1.000
84.0 53 1 0.666 0.16005 0.416 1.000
89.0 52 1 0.657 0.16326 0.404 1.000
95.0 51 1 0.648 0.16648 0.392 1.000
99.0 50 1 0.639 0.16972 0.380 1.000
101.0 49 1 0.630 0.17293 0.368 1.000
109.0 47 1 0.621 0.17611 0.356 1.000
148.0 45 1 0.611 0.17927 0.344 1.000
152.0 44 1 0.601 0.18236 0.332 1.000
164.0 43 1 0.592 0.18551 0.320 1.000
185.0 41 1 0.583 0.12737 0.380 0.894
187.0 40 1 0.574 0.12889 0.370 0.891
206.0 39 1 0.565 0.13036 0.359 0.888
218.0 38 1 0.556 0.13180 0.349 0.885
262.0 37 1 0.546 0.13320 0.339 0.881
284.0 35 2 0.527 0.13585 0.318 0.874
307.0 33 1 0.517 0.13707 0.308 0.869
333.0 32 1 0.507 0.13823 0.297 0.865
339.0 31 1 0.497 0.13930 0.287 0.861
342.0 29 1 0.486 0.14029 0.276 0.856
583.0 21 1 0.471 0.14187 0.261 0.850
674.0 17 1 0.452 0.14361 0.243 0.843
732.0 16 1 0.433 0.14506 0.225 0.835
851.0 14 1 0.410 0.14622 0.204 0.825
979.0 11 1 0.383 0.14698 0.180 0.813
995.0 10 1 0.356 0.14735 0.158 0.801
1031.0 9 1 0.330 0.14743 0.137 0.792
>
> # These should all give the same answer
> # When there are offsets, the default curve is always for someone with
> # the mean offset.
> j.age <- jasa$age -48
> fit1 <- coxph(Surv(futime, fustat) ~ j.age, data=jasa)
> fit2 <- coxph(Surv(futime, fustat) ~ j.age, jasa, init=fit1$coef, iter=0)
> fit3 <- coxph(Surv(start, stop, event) ~ age, jasa1)
> fit4 <- coxph(Surv(start, stop, event) ~ offset(age*fit1$coef), jasa1)
>
> s1 <- survfit(fit1, list(j.age=fit3$means), censor=FALSE)
> s2 <- survfit(fit2, list(j.age=fit3$means), censor=FALSE)
> s3 <- survfit(fit3, censor=FALSE)
> s4 <- survfit(fit4, censor=FALSE)
>
> all.equal(s1$surv, s2$surv)
[1] TRUE
> all.equal(s1$surv, s3$surv)
[1] TRUE
> all.equal(s1$surv, s4$surv)
[1] TRUE
>
> # Still the same answer, fit multiple strata at once
> # Strata 1 has independent coefs of strata 2, so putting in
> # the other data should not affect it
> ll <- nrow(jasa1)
> ss <- rep(0:1, c(ll,ll))
> tdata <- with(jasa1, data.frame(start=rep(start,2), stop=rep(stop,2),
+ event=rep(event,2), ss=ss, age=rep(age,2),
+ age2 = (rep(age,2))^2 * ss))
> fit <- coxph(Surv(start, stop, event) ~ age*strata(ss) + age2, tdata)
> # Above replaced these 2 lines, which kill Splus5 as of 8/98
> # Something with data frames, I expect.
> #fit <- coxph(Surv(rep(start,2), rep(stop,2), rep(event,2)) ~
> # rep(age,2)*strata(ss) + I(rep(age,2)^2*ss) )
> all.equal(fit$coef[1], fit3$coef)
[1] TRUE
> s5 <- survfit(fit, data.frame(age=fit3$means, age2=0, ss=0), censor=FALSE)
> all.equal(s5$surv[1:(s5$strata[1])], s3$surv)
[1] TRUE
>
>
>
> proc.time()
user system elapsed
1.011 0.101 1.104